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ABSTRACT 


Title of Dissertation: Stability and Interaction of Coherent 

Structure in Supersonic Reactive Wakes 

Suresh Menon, Doctor of Philosophy, 1983 

Dissertation directed by: Dr. John D. Anderson, Jr. 

Professor 

Department of Aerospace Engineering 
and 

Dr. Shih I Pai 
Professor Emeritus 
Institute of Physical Science 
and Technology 

A theoretical formulation and analysis Is presented for a study of 
the stability and interaction of coherent structure in reacting free shear 
layer. The physical problem under investigation is a premixed hydrogen- 
oxygen reacting shear layer in the wake of a thin flat plate. The coher- 
ent structure is modeled as a periodic disturbance and its stability Is 
determined by the application of linearized hydrodynamic stability theory 
which results in a generalized eigenvalue problem for reactive flows. De- 
tailed stability analysis of the reactive wake for neutral, symmetrical 
and antisyimietrical disturbance is presented. Reactive stability criteria 
is shown to be quite different from classical non-reactive stability. The 
interaction between the mean flow, coherent structure and fine-scale tur- 
bulence is theoretically formulated using von-Karman integral technique. 
Both time-averaging and conditional phase averaging are necessary to sep- 
arate the three types of motion. The resulting integro-differential 
equations can then be solved subject to initial conditions with appropri- 
ate shape functions. In the laminar flow transition region of interest, 
the spatial interaction between the mean motion and coherent structure is 



calculated for both non-reactive and reactive conditions and compared with 
experimental data wherever available. The fine~scale turbulent motion is 
determined by the application of integral analysis to the fluctuation 
equations. Since at present this turbulence model is still untested, 
turbulence is modeled in the interaction problem by a simple algebraic 
eddy viscosity model. The applicability of the integral turbulence model 
formulated here is studied parametrically by integrating these equations 
for the simple case of self-similar mean motion with assumed shape func- 
tions. The effect of the motion of the coherent structure is studied 
and very good agreement is obtained with previous experimental and theo- 
retical works for non-reactive flow. For the reactive case, lack of 
experimental data made direct comparison difficult. It was determined 
that the growth rate of the disturbance amplitude is lower for reactive 
case. The results indicate that the reactive flow stability is in quali- 
tative agreement with experimental observation. 
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I. INTRODUCTION 


The study of the interaction between turbulence and chemical reac- 
tions is a region of continued interest in fluid dynamics. A better 
understanding of the complexity involved will enhance design and predic- 
tive capability of advanced and more efficient propulsion systems. Of 
particular interest here is the combustion system of a Supersonic Combus- 
tion RAM Jet engine (SCRAMJET) being currently studied at NASA Langley 
Research Center [1]. In the combustor of a SCRAMJET, hydrogen fuel is 
injected at near sonic condition into a supersonic air stream and com- 
bustion occurs in the recirculatory zone inside the combustor. The analy 
sis of such a high temperature, turbulent reactive flow is complicated by 
the ongoing interaction between fluid mechanical and chemical effects, 
whereby, the reactants first mix and then combine to release chemical 
energy which in turn significantly alters the flowfield. The analysis 
of such an interaction is so complicated that without some simplifica- 
tion the problem remains intractable. However, with some reasonable 
assumptions an understanding of the complex interaction can be obtained. 

The application of the eddy-viscosity models developed for non- 
reactive flows in combustion studies has met with only limited success 
[2]. The use of gradient diffusion models have been shown to be incorrect 
in reactive and recirculatory flows where counter gradient diffusion is 
present [3, 4]. It became quite clear that new approaches were necessary 
to handle the complexity of reactive flows. This led to the counter gra- 
dientmodel [3] and probabi 1 ity distribution and functional formulations [5- 
8]. Though these models have been successful in predicting some of the phe 
nomena inherent in reactive flows, there are still many unanswered ques- 
tions concerning the interaction between chemical kinetics and flow 
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turbulence. 

Ever since the existence of large scale coherent structures In tur- 
bulent shear flows became evident, the effect these coherent structures may 
have on combustion processes has been subject to speculations [9-12]. This 
has led to various new models, for example. Spalding's ESCIMOmodel [13,14] 
which accentuates a genuine property of turbulent flows that is evidently en- 
hanced by large coherent structures: the stretching of local flame elements 
due to the straining motion of the flow. Anothermodel by Marble and 
Broadwell [15] studies the diffusion flame structures as opposed to the 
premixed flame model of Spalding. The Influence of fluid dynamics on com- 
bustion has been considered by Chorin [16] using a numerical method based 
on vortex dynamics which has been also extended by Ghonlem, Chorin and 
Oppenhelm to non-constant density flows with the aim of modeling combus- 
tion in coherent structures [17]. 

There Is now much experimental evidence of the existence of large scale 
coherent structures In free shear turbulent flows [9,10,18]. Through flow 
visualazatlon, Moore [19] showed that a turbulent round jet also has a defi- 
nite coherent structure that starts as an Instability wave In the shear 
layer. Earlier experiments by Pal [20] had first pointed out the existence 
of secondary flow Inside rotating cylinders. More recently, Gan j1 and Sawyer 
[ 21 ] observed large structures dominating mixing layer that develops behind 
a step under non-reacting and reacting conditions. In comparing reacting 
and non-reacting flows, they found that the reacting eddies have a lower 
growth rate, and more closely distributed In space and have a slightly 
smaller ratio of coalescence than non-reacting eddies. In turbulent 
flames Vule et. al [22,23] found that combustion driven Instabilities 
effect the coherent structure growth and decay. In fact, they found 
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at least two combustion driven instabilities, an Inner high frequency 
and an outer low frequency phenomena. These instabilities are genuine 
properties of flames and do not occur in non-reactive flows, leading to 
the term 'combustion driven coherent structures' [24]. Such experimental evi- 
dence supports Roshko' s [25] conclusion that coherent structures play a cen- 
tral role in the development of many turbulent shear flows such as mix- 
ing layers, boundary layers, and the early regions of jets and wakes. 

It has also been noted that combustion in non-premixed flames seems 
to conserve coherent structures in the flow by delaying transition to 
turbulence [26]. All available data seems to indicate that coherent struc- 
tures are potentially more important in combustion systems than any other 
flow systems due to the strong influence they have on the turbulent mix- 
ing of reactants and to the stabilization of existing structures by com- 
bustion. 

The transport processes across the mixing layer is considerably 
under-predicted by all theoretical models leading to the point of view 
that the effect of large scale coherent structures on scaler mixing pro- 
cesses cannot be predicted by methods that use scalar flux approximations. 
Therefore, it seems clear that a different closure model is necessary to 
handle the combustion problem in turbulent flow. 

The present investigation considers the theoretical analysis of the 
laminar-turbulent transition of compressible reactive wakes. Experimen- 
tal measurements in non-reactive wakes behind flat plates and slender 
wedges [27,28] have shown remarkable similarity with the low-speed wake 
transition analysis by Sato and Kuriki [29]. Due to their inherent dynamic 
instability, wakes sustain travelling wave disturbances. The development 
of these instability waves and their consequent interaction with the mean 
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velocity, thermal and concentration fields and the fine-scale turbulent 
fields constitutes the interaction problem considered here. The coher- 
ent structure discussed above is modelled as an instability wave which 
developes from a linear growth region into a nonlinear growth and finally 
into three dimensionality. The disturbance amplitude is very small in 
the linear region and the mean field is uncoupled from the disturbance 
field. However, in the nonlinear region, the amplitude becomes large 
and there is a strong interaction between the mean field and disturbance 
field causing the mean field to decay more rapidly than in the linear 
region. Beyond the nonlinear region, the disturbance becomes three 
dimensional and for high enough Reynolds number the flow becomes turbu- 
lent. This is the general picture of wake transition although the actual 
extent of each region depends upon flow field parameters like the Rey- 
nolds number and Mach number. 

This analysis considers the motion of the turbulent fluid as a com- 
bination of three distinct motions: the mean motion, the large-scale 
coherent structure motion and the fine-scale turbulence. Such a splitting 
procedure was first used by Reynolds and Hussain [30] and has been used 
extensively by Liu et. al [31-34] in their study of coherent structures. 
More specifically, Liu and Merkine [33], Liu and Alper [34], Alper and Liu 
[35] and Gatski and Liu [36] have studied the interaction between a mono- 
chromatic component of the large-scale coherent structure and the fine- 
grained turbulence in developing mean flows with inflexional profiles. 
There, the nonlinear interactions between the three components of flow are 
depicted in terms of the non-equilibrium adjustments between the mean shear 
layer growth rate and the integrated energy densities of the large-scale 
structure and the fine-grained turbulence. Their analysis was limited to 
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Incompressible non- reactive flow, however, they found reasonable agree- 
ment with experiments for the spatial spreading rate of the shear layer. 

The present theoretical study is an analysis of the nonlinear inter- 
action between the three motions in a compressible multicomponent reac- 
tive flow. Two averaging procedures are necessary to consider the inter- 
action among the three components of flow: conditional phase averaging 
which explicitly filters the coherent structure from the total fluctua- 
tions containing both random and coherent components; and the conven- 
tional Reynolds averaging which separates the mean flow from the fluc- 
tuations. Due to the complexity of the resultant equations, finite- 
difference computations of the interaction is very complicated. Stuart 
[37] and Ko, Kubota and Lees [38] analysed the development of finite 
amplitude disturbances in incompressible shear flows using von Karman 
integral formulation. Integral considerations were also used by Liu 
and Gururaj [39] to study compressible wake transition and they obtained 
good agreement with experimental data for hypersonic wakes [27, 40]. Ac- 
cordingly, the present study incorporates von Karman integral method to 
obtain the conservation equations for the mean flow. The mean flow is 
then characterized by the wake width, the wake centerline values of the 
mean velocity defect, the mean temperature excess and the mean specie 
mass fraction defects. The disturbance is characterized by its am- 
plitude and its variation across the shear layer is determined by the 
application of hydrodynamic stability theory to the disturbance equa- 
tions. The integral technique is also applied to the fine-scale turbu- 
lence equations to achieve closure. Though the present formulation is 
simplified by various assumptions (to be discussed later), it is expected 
that this formulation can be used to obtain a better understanding of the 
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Interaction between the three components of flow during combustion, an 
area of research that has not been studied so far. 

The physical problem studied at present is the growth of a multi - 
component free shear layer in the wake of a flat plate (Figure 1). For sim- 
plicity we consider a supersonic premixed stream of hydrogen, oxygen, 
hydroxyl radical and water vapor in the shear layer. The interaction 
problem for both non-reactive and reactive flow conditions is studied 
and the changes in the organized motion due to reaction is discussed. 

In Chapter II the governing conservation equations of the interaction 
problem is formulated and the various assumptions used are discussed. 
Chapter III describes the numerical methods used to solve reactive sta- 
bility problem and the integral equations of motion. In Chapter IV the 
results of the numerical calculations are presented and compared with 
other theoretical studies and experimental data wherever possible. 









II. FORMULATION OF THE PROBLEM 
2.1 The Governing Equation 

The general conservation equations for a multicomponent mixture in a 
cartesian coordinate system can be written in dimensional form as 
a) Mass conservation 


3p + — i _ 0 


3pu. 


( 1 ) 


b) Species Conservation 


3pY. 


~9l ■ 8x. ^'k 

c) Momentum Conservation 
3 




k = 1 • • *v 


3pu. ^ 3 
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d) Energy Conservation 
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(3) 


|^(pH-p)+^(pu.H-u.T.^+q^) = 
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e) Equation of State 

p = pTJR|^Y|^ 
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Here t.. is the shear stress tensor given as 
^ 0 


■'id 


-TP 


«<< * P <73^^ 857' 


3 ^ ^ “ij 
H is the total enthalpy defined as 


H ' 1 (u^ + v^) + Z (n^T + 
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-)R. gives the translational and rotational de- 
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grees of freedom (by njj) and 6* gives the i-th vibrational degree of 
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freedom for the k-th specie and ej* is the heat of formation for k-th 
specie. 

is the heat flux defined as 

( 8 ) 

where h^ = e^^ + T Is the static enthalpy for the k-th specie; ej^ and 

Rj^ are the specific internal energy and specific gas constant for the 

k-th specie, u. is the diffusion flux velocity for the k-th species 
^i 

given by 

°k '^k 

v = -r,^ ’ ^ = (9) 

where is the molecular diffusivity of the k-th specie and k* is the 
thermal conductivity. 

2.2 Filtering Procedure 

The above equations (1 - 9) are now reduced to the integral form to 
be used in this analysis. Though turbulence is essentially three-dimen- 
sional, for numerical simplicity, we consider two-dimensional flow. Ex- 
tention to three-dimensional flow will not change the governing equations 
appreciably. In wake type flows, it has been shown that the detailed 
distribution of the disturbances is not very sensitive to the viscous 
terms leading to the 'inviscid' considerations [39,41]. Though molecular 
transport phenomena is important in turbulent reactive flows, it is 
smaller in comparison to the transport due to turbulent stresses and may 
be neglected in comparison. We include only a simplified form of the 
dissipation terms in the governing equations and neglect the fluctuations 
in the transport properties. 

Any instantaneous flow variable, q(x,y,t) is then decomposed as 
follows: 
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q(x,y,t) = q{x,yj + q'(x,y,t) + q"(x,y,t) 


no 


Here q is the mean flow component, q' and q" the corresponding 

coherent structure and fine-scale turbulent components respectively. 

On using (10) in the Navier-Stokes equations and time averaging we get 
the equations governing the mean motion. Subtracting the mean motion 
equations from the original equations and then taking phase averaging 
we extract the equations governing organized motion. When we subtract 
the equations for organized motion from the total fluctuation equation we 
obtain the equations governing fine-grained turbulence (see Appendix B). 

The two averaging processes as discussed above are defined as follows; 

Time average is 


whe re T 
at least. 


q(x,y,t)dt 


( 11 ) 


0 

is greater than the period t of the large-scale structure 
The conditional phase average is 


<q(x,y,t)> = 


Lim X 

f^co N 


N 

t q(x,y,t + m) 
n=0 


( 12 ) 


The conditional and time average of a turbulent quantity, {y. ,y ,t) , 

are zero by definition. The conditional average of a large scale structure 
quantity, q'(x,y,t) reproduces itself and its time average is zero. We 
assume that the two components of the fluctuations are not correlated. 
Furthermore, the conditional average of two fine-scale turbulent quantities 
after subtracting the steady part, <q"q"> - q"q" , is periodic and 
oscillates at the same frequency as the large-scale structure [31]. 

Though the large-structure may have many frequency components, at 
present we consider the propagation of only the fundamental component. 
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This Is a reasonable approximation because it has been shown experiment- 
ally that through proper control it is possible to get a single mode 
propagation in the flow [42], Furthermore, as a first approximation, the 
energy this fundamental mode exchanges with other frequency components 
can be neglected with respect to the energy it exchanges with the mean 
flow or the fine grained turbulence. 

2.3 Conservation Equations for the Mean Flow ‘ 

To derive the mean flow equation we nondimensionalize all variables 
with respect to the free stream conditions which are assumed to be con- 
stant. Thus, velocity and coordinates are nondimensional ized by the free 
stream velocity, u , and reference length, L respectively; the pressure, 
temperature and density are made dimensionless by their corresponding 
free stream values p^, T^ and p^. 

The time averaged equations obtained are complicated and contain 
many unknown correlations of the fluctuations. To reduce the complexity, 
we make some approximations which have been shown to be reasonable by 
past analysis of wake flows. Since we are considering free shear 
flows at present, we can apply the boundary layer approximations without 
losing significant accuracy. We further use von Karman integral technique 
to integrate the governing equations across the shear layer (normal direc- 
tion) and obtain a set of integro-differential equations. In deriving 
these equations we assume that all fluctuation correlations vanish far 
away from the flat plate and shear layer region. The problem then reduces 
to determining the shape functions for the mean flow variables associated 
with the von Karman integral formulation. 

After some manipulation we get the following integral equations: 
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2.3.1 


Specie Integral Equation 


^ f P p(Vr’) 


+ R*; 
dx k 


k = 1 


V species 


(13) 


Here, ^ and IT are the dimensionless mean density and streamwise 
velocity respectively. Yj^ is the normalized mean mass fraction of the 
kth specie defined as 

fraction and Y, is the freestream value of mass fraction of the kth 
k® 

specie . 


Ip, is the species diffusion flux integral for the kth specie 

^ c 

I ^ 

made up of contributions from the organized structure motion, D 

I ^ 

and the fine- scale turbulence, . They are given as 


[p~iT (Y|^ - 1 ) + p u' Y|^' + u p'Yj^' ] dy 


(14) 


and 


so that 



= I* [p"u" (Y^-1) + p u"Y^' + u p"Y^" + p"u"Y^’] dy 



(15) 

(15) 


is the integrated mean rate of production of the kth specie given 
as < 

R* . fir- dy on 

»CB 

Here 71^** is the dimensionless mean rate of production of the kth specie 
and is defined in Appendix B. 
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2.3.2 Streamwise Momentum Integral Equation 


dx 


p u(u - 1 ) dy = - 


dl 


M 


dx 


The momentum diffusion integral, is defined as 

'm ' ’h' * " 'h 


where 


= I [p ' (2u-l) + p u'^] dy 

— OD 

= j 17 ^ (2U - 1) + ^ u"^ + p"u"2] dy 



( 18 ) 


(19) 


(20) 

( 21 ) 

( 22 ) 


dl 


M. 


dx 


In the free shear layer, pressure is nearly constant and so the term 


is negligible in equation (18) and 


2.3.3 Normal Momentum Integral Equation 


P = 1 - ^ ^ 


{p (u'v' + u"v'*) 


+ u (pV + + pW^ } dy] 


(23) 


It has been shown by Liu and Gururaj [39] that for wake flows the last 
term on the right side of equation (23) is negligible relative to the 
others and pressure p can be determined from the reduced equation. 
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2.3.4 Energy Integral Equation 


d 

'3^ 


r 

— OD 


P u ( H - 1 ) dy = 



( 24 ) 


The total enthalpy integral, is defined as 

If = Ij" * (25) 

where 


^00 

= I + u 7 h^+ (H- 1 ) TTH dy (26) 

— 00 

T f“ 

= J [p u H^' + u p"H" + (H-1) p"u" + p"u"H"] dy (27) 

— OP 

and the mean total enthalpy, iT is given as ' 

H = + (28) 


2 

where Q 
energy. 


= ^ (u^ + + 

The fluctuation 


v' + u"^ + v"^) is the total mean kinetic 
contribution to the total enthalpy are 




+ E 
k=l 



(29) 


H" . uu". 1 . T“7^ . L„^) 


+ E 
k=l 



(30) 
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The term R and L, , appearing in equations (29) and (30) are defined 
explicitly for any variable (j) and \p as 


+«!•>"> - <J'>" 
and 


C, , D. and e. are constants given as 
1 K 


CT V, 


<^1 ' 


P U 
^ 00 00 


k - SITT 




kf K„ 


(31) 


(32) 


Here, R is the universal gas constant and ^ jand. are the molecular 

weight and heat of formation of the k-th 'i^eci^ /' Also 

N 0 

N./ I 

-] 


'k./T 

\ ^ ^ -?T“ 



(33) 


i=l (e ’^i''^-!) 


gives the translational, rotational (by Nj^) and vibrational degrees of 

freedom with 6, as the characteristic dimensionless temperature for 
^i 

the i-th vibrational frequency of k-species. 


2.3.5 Equation of State 


p = E C* [p T Y^+ F(TT^ + T"Y\) + T(^V^ + 
k=l 


(34) 


+ ^P'T' + 7T^) + p"T"r^J 


Where 


p RT 

r* = — Y 

k M|^p„ k„ 


(35) 


2.3.6 Mean Kinetic Energy Integral 

It is also necessary to derive the integral equation governing the 
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total mean kinetic energy and can be written as 




d I 


— ^ - I + I - 
dx ^RS <*> 


(36) 


Here, represents the mean flow kinetic energy diffusion flux integral 

and is defined as 


^kE " ^kE ^kE 


(37) 


c T 

where and are the diffusion flux integrals due to organized 

structure motion and fine-scale turbulence respectively. They are given 
as 


IkE^^ = (3u^-l) p'u' + p u u'^] dy 
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(38) 


and 


^kE 


l\ (3u^- 1) p"u" + p u u"^ + p"u"^] dy 


(39) 


Ic_ represents the total turbulent stress production integral contributing 
R5 

to the mean flow. It is defined as 


c T 
^RS " ^RS ^RS 


(40) 


where is the Reynold stress production due to organized structure 

Kb 


motion 
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^RS 
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and the corresponding stress production integral due to fine scale 
turbulence 


T ^ ^ 

^RS = - J ^ + P"u"^)|^ + {p TV" 

-00 

- V Tv" + pV'v") + (p V"^ + p-V'^) 1^] dy 

Ip denotes the pressure integral given as 

I = I + I 
P P] P2 


(42) 


(43) 


where 


00 



I = -C 
p^ ^1 dx 

00 

j (p - 1) u dy 

(44) 

and 

_co 



'P2 = f 

(p.,) (3U . 3l) ,, 

(45) 


Here, I is the pressure diffusion integral and I gives the work done 
Pi P2 

due to expansion. The mean flow viscous dissipation integral is simply 

taken to be 

j ^ 0 ' <«> 

—00 

v/here p is the dimensionless viscosity coefficient and Re, = p u L/p 

00 00 00 

is the Reynolds number based on freestream conditions. 


2.4 Shape Assumptions for Mean Flow 

Experimental investigation indicates that in a wake the momentum 
thickness 6*, is very nearly constant. Equation (18) may then be rewrit- 
ten as 
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(47) 


00 

/ p u (u - 1) dy = Constant 

-CD 

by neglecting the integral . 

On using the definition of momentum thickness 

e* * p u (1 - u) dy = T- = constant (^8) 

0 ^ 

where, Cp is the body drag coefficient, along with equation (47) a rela- 
tion between the wake width and the centerline velocity defect can be 
obtained as shown later. 

The solution of the governing integral equations for the mean motion 
depends upon the proper choice of the shape functions characterizing the 
variation of the mean flowfield in the wake. It is well known that the 
Von Karman integral method allows quite accurate calculations of the 
gross features while remaining quite insensitive to the finer-details, 
provided of course, that the shape function representation is physically 
reasonable. The distribution of the mean flow variables across the wake 
has been experimentally shown to be nearly Gaussian [27.65]. Therefore, the 

the mean flow is represented by 

u (x.n) = 1 ■ Vj. (x) exp (-nM 

T (x.n) = 1 + ■^c 

Y. (x,n) = 1 (x) exp (-n^). a = l.---v 

where all variables are non-dimens 1 onali 2 cd with respect to their edge 
values. Here, V (v), T (x) and (x) are the mean horizontal velocity 
defect, temperature and k-specie mass fraction excesses at the centerline 
of the wake. The coordinate n, is defined as the Howarth transformed 
variable given by 


17 



( 50 ) 


y _ 

n = / p dy / 6 (x) 

where 6(x) is the dimensionless shear layer thickness which can be shown 
to be related to the centerline velocity defect V^(x) by equation ( 48) 
to give 

' V^(Ci CzV^) 

« 2 „ 2 
where Ci = / e ^ dn and, C? = /“e^'^ dn 

»oo ^ ■•oo 

Therefore, the mean flow can then be characterized by V (x), T (x), 

w V» 

and Y. (x). The variation of mean density p(x, n) and pressure, 

_ c 

p(x, n) can be represented in terms of V^(x), (x) and T^(x) in the 

c 

governing equations with the help of normal momentum equation (23) and 
the equation of state (34). The mean flow field can then be completely 
defined in terms of the variable V^(x), T^(x) and Yj^ (x). 

2.5 Coherent Structure Closure 

The equations described in Section 2.3 govern the mean flow. How- 
ever, they contain correlations due to the fluctuations of the coherent'-' 
structure variables and also fine-grained turbulence. Here, we discuss 
the closure of the coherent structure correlations. As mentioned before, 
the equations governing the coherent structure motion are derived using 
the hydrodynamic stability theory. The coherent structure is assumed to 
have a wave-like periodicity as is defined by its amplitude A(x) and 
its frequency of oscillation 8 . By choosing appropriate shape functions 
for the coherent structure variables we derive the local eigenvalue equa- 
tion which is solved using local linear stability analysis. 

To facilitate numerical solution we define local coordinates given 
as 
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x-x^ 


( 52 ) 


where for the eigenvalue problem the dimensionless shear layer thickness 
6(x) is considered locally a constant and is a location close to x 
such that the mean flow streamwise gradients are negligible. 

2.5.1 Shape functions for Organized Structures 

The shape functions assumed for the coherent structure in their 
dimensionless form can be written in the local coordinates as 


u'(C,n,T)/V^ 


2{ti) 



av(n) 
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p(n) 

, p’ (C.n »t) 

= A(x) 

P(n) 
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T{n) 

Y|^(C.n,T),k=l,v 
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-<> 

II 


c . c + 0([A|^) 


(53) 


Here, c 'C denotes the complex conjugate and the complex eigenfunctions 
of the fundamental component u, v,p,p,T and Y|^ correspond to the 
wave streamwise velocity, normal velocity, density, pressure, temperature 
and kth species mass fraction, respectively. Furthermore, p = p* 

CD 

is the real dimensionless frequency, p* the physical frequency; and 
o = a* 6(x)L is the dimensionless complex wavenumber and a* its dimen- 
sional value. The local stremwise derivatives of the disturbances 
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are evaluated by the local relation [31]. 

^=iaA (54) 

On using (53) and (54) in the linearized inviscid equations for the 
organized motion and on some algebraic manipulation (see Appendix A), the 
governing equation in terms of the pressure perturbation p(n) can be 
written as 


P" 


+ p- [. 


2w’ 

(w-c) 




- (w-c)\2 ^ 

+ p [ V 

T 


A 1 1 ^ T 1 2 c 2^2 _ Cl 

^C-|S " a S T - 0 


(55) 


where w = (u-l)/V^ and c is the dimensionless complex phase velocity 
defined as c=(3/a-1 )/V^. Here, prime denotes differentiation with re- 
spect to n» and 


S = z'’ C* , S' = z'’ C* and M = z'’ m^D^T^ (56) 

k^l k“l k~l 

Also, v(= Vp + ivj) and "3(= Ip + ilj) are complex terms appearing 
due to the presence of finite rate kinetics and are described in Appendix 
A. We note here that this eigenvalue problem eq. (55) reduces to the 
classical Lees and Lin [43] problem when the flow is non-reactive and 
single specie. For multi -component non-reactive flow 


* 

Vn 


= 1 . 




(57) 


The eigenvalue equation (55) is solved for the assumed mean flow 
profiles and appropriate boundary conditions for a given frequency 6 
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and complex wave number a. All other eigenfunctions u, v, p, T and 

^ A 

Y|^ are subsequently determined from the local solution for p(n)- 
Inherent in such a consideration described above is the assumptions of 
'local similarity' which Implies that the eigenfunction p(n) adjusts 
instantaneously to the mean flow. To complete the closure of the 
coherent structure we need to determine the amplitude A(x) appearing 
in equation (53). This is given by the equation governing the time- 
-averaged kinetic energy production associated with wavelike motion of 
the organized structure. 


2.5.2 Coherent Structure Kinetic Energy Integral 
c 1 ^ ? 

where k = (u' + v* ) is the organized structure mean kinetic energy. 

2 

The eigenfunctions are normalized such that lA| is the dimensionless 
kinetic energy per unit length, of the instability wave across a slice 
of the shear flow [39]: 


00 

W _oo 


k dn = 
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00 

I (u'^ + v'^) dn 


( 59 ) 


C -00 

and may be considered as the energy density in terms of the transformed 
coordinate n. 

Ij^^ is the Reynolds stress production integral defined in equation 
cT 

(A1). Ij^j is the kinetic energy integral for the exchange between 
the organized structure and the fine scale turbulence and is given as 
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is the pressure gradient work integral of the disturbance and is 


given as 
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and the viscous dissipation integral for the disturbance is simply 
taken to be [41]: 
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I,^ generates mean flow thermal energy at the expense 
kinetic energy. 


of disturbance 


The energy exchange integral which describes the interaction 

between the coherent structure and fine-scale turbulence is in a form 
that requires closure. The terms R^^^ are specified in terms of an eddy 
viscosity and disturbance strain rate model similar to Liu [31] such that 
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On neglecting the terms larger than 0 (|a!^), can be written as 
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The eddy viscosity f,€ used in this closure is obtained in terms of 
the turbulent kinetic energy and dissipation determined in the next sec- 
tion. and Op* is a constant. To first approximation, ^ • 

2.6 Closure for the Fine-Scale Turbulence 

The equations governing the fine-scale turbulent quantities can be 
obtained by subtracting the contributions to the mean and organized 
structure motion from the instantaneous equations. These equations 
are quite complex and will not be shown here. Here we present a closure 
model for the turbulent fluctuation correlations. Since the mean motion 
and the organized structure motion have been formulated using integral 
technique, we follow similar steps to develop the turbulent equations. 
Essentially, we manipulate the turbulent fluctuation equations and write 
equations for turbulent-density correlation . turbulent-specie mass 
fraction correlation turbulent - total enthalpy correlation H" 
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and turbulent kinetic energy q" . On integrating the resultant equations 
across the shear layer and assuming that all fluctuation correlations 
vanish far away from the shear layer we get the following integral 
equations. 

2.6.1 Density Fli^tuation Integral Equation 
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I dx 
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dy = - I + I - I 
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where I is the production of p" due to mean strain rate and is given 


by 
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Ip is the diffusion flux integral given as 
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and Ip is the! di ffusion integral for higher order correlation. 
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We observe that there are terms given above that require closure assump- 
tions. Hence we make the following assumptions. The turbulent density 
strain rate is given as 

p" (^ -^ ^) = ^ 

where e is a constant and \ is the Taylor's microscale characteristic 

p 

of this problem to be defined later. 
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We further use a simplified gradient formulation to model terms 

M 


likeL,. Thus 



(72) 
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The integrals I and I in equations (53) and (54) then become 
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where S and S are constants that have to be adjusted with the help 
p p-| 

of experimental data and is the eddy mass diffusivity. 

2.6.2 Specie Fluctuation Integral Equation 


where I is the convective flux integral given as 
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Ip is the diffusion flux integral due to mean concentration gradient, 
*^k 
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( 78 ) 






I is the turbulent diffusion of species and is given by 
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kj dx 
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I is the production integral for the k-th species due to chemical 
Pk 

reaction „ 

’p - [ ‘®°’ 

^ wOO 

where rj,®" is the turbulent production rate of the k-th specie due to 

finite-rate kinetics and is obtained from the total species production 

rate. Ip is the diffusion flux integral due to species diffusivity. 

• „ 

[ P D, Y," dy (81) 

vOD 

Here, is the mean molecular diffusivity of the k-th specie and 
is defined later. We neglect any fluctuations in Dj^ at present, 

To model the unknown correlations we follow similar approximation 
steps. For example. 
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On replacing these approximations in the integrals in Equations (61), 
(62), (63) and (64) we get 
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where S|^ and Cq are constants that are adjusted with available experi- 
mental data. 
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2.6.3 Turbulent Kinetic Energy Integral Equation 
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where 1 q"^ = j (u"^ + v"^) is the turbulent kinetic energy and 

is the Reynold stress integral associated with the fine-scale turbulence 

interaction with mean flow and is given in Equation (33). 

I is the Reynold stress production integral associated with 
RS 

the interaction of fine scale turbulence with itself and is given as 
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On neglecting fourth order terms at present and following the previous 
closure approximation we get 
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where k^* and k^* are adjustable constants, 
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The terms appearing In (90 ) measure the effect of viscous decay on the 
second order correlations. We follow here Varma et. al. [44] approach 
which is an anistropic dissipation model and is shown to have correct 
expression in the limits of small and large Reynolds number. In general 
we write 




(91) 


where X is the Taylors microscale length and is related to the macro- 
scale A by the Rotta's form 






(92) 


where 


A as 0.5 6 and v = v /'p 

The study of incompressible boundary layer led to values of constants 
a and b as a = 3.25 and b = 0.125. However, for compressible flows 
they may have to be readjusted [44]. 

The dissipation correlations involving scalar fluctuations are 
modeled in the following form 
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Ip is the pressure work integral due to fine scale turbulence 
and can be split into two parts. 

P P] P 2 


( 94 ) 


where I is the pressure work due to expansion and is given as 
Pi 
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and I is the pressure-strain rate 
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To determine these terms we follow Varma et. a1 [44] approach to model the 
pressure diffusion term as follows 

(97) 




p"u" = P 


A IP, (q"')’/' * P, A 1^] 


where p^ and are constants. We take P 2 = 0*1 and p^ which is 
associated with the mean strain, needs to be adjusted with available 
data. The pressure-strain rate is given by the form developed by Launder 
et. al . [45] which includes the tendency towards isotropy terms. 
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where the first-term is proportional to the anistropy of the turbulence, 
The other terms are due to the interaction between turbulence and mean 
flow and is given by production terms 
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Onrp thp di??iDation 


depends upon the size of eddy we apply the knowledge from simple spectral 

analysis to find the variation of dissipation rate in the shear layer. 

We shall return to the determination of e later on. 

Finally, is the turbulent viscous dissipation integral due 
♦ 

to fine scale turbulence. We model this term in a form similar to 
the viscous dissipation due to organized structure motion. It is then 


given as 
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where y is the mean molecular viscosity and we neglect any fluctuation 
in transport properties. generates mean flow thermal energy at 

the expense offine-scale turbulent kinetic energy. 

Since the energy conservation requires balance of the kinetic ener- 
gy in the three modes of motion, Mankbadi and Liu [48] determined equations 
where the stress term ii^ equation (58) appeared with an opposite 

sign in equation (87). Since at present, equation (87) does not con- 
tain any term defining the interaction between the turbulent and dis- 

. . 1 T TC . T 

turbance field, to the first approximation we replace 
in equation (87). 

2.6.2 Total Enthalpy Fluctuation Integral Equation 
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where lu is the convective flux integral due to mean strain rate and 
is given by 
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(107) 
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Ijj Is the diffusive flux integral due to mean total enthalpy gradients. 
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and is the turbulent ehthalpy diff usion in tegral 
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Im is the dissipation integral which we neglect at present. 
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Following the same closure approximations we get the following 


integrals 

'h. 


OD 

-- [[piTTMu-Su ^)?TF]|?dy 


X 


3X 


Ob 

■j 


[p dy 


and 


■-li ] f 


- 1^) dy 


(110) 


(111) 


where S„ and Sy are constants. 
H] 
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Here, ^ Is the turbulent total enthalpy diffusivity and Is defined by 
the relation 

“ Si ^ 

where Cu is a constant and e„ is the dissipation rate for total 

-? 

enthalpy fluctuation correlation H" . , 

It should be noted here that the closure form used here in terms of 
contains implicitly the combined effect of turbulent kinetic energy, 
the turbulent thermal energy and the turbulent chemical energy. There- 
fore, the turbulent total enthalpy diffusivity „ can be used instead of 
D,, in the equation for Y. . Since, momentum dissipation e and enthalpy 
dissipation are very important in the turbulent field development, we 
determine them by using a simplified spectral analysis using the multiple 
scale approach of Hanjalic et al [49]. The underlying physical mechanism 
of energy transfer in turbulent flows is complicated further due to the 
presence of finite rate kinetics and therefore a multiple scale process 
of energy transfer seem more appropriate for the present problem. 

2.6.5 Basis for Spectral Analysis 

It is well known that turbulent energy transport occurs due to the 
motion of turbulent eddies in the flow. Besides kinetic energy of the 
motion, thermal energy and chemical energy due to reactions is also 
being transported when different species mix and react at the molecular 
level within these eddies. 

The complex interaction of the turbulence and chemical reaction 
leads to growth and decay of turbulent eddies in a manner significantly 
different from non reactive flows. Previous analysis of energy cascade 
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phenomena In turbulent flows assumes that energy is extracted from the 
mean flow and transfered in proportion, from the largest eddies to the 
smallest eddies where it is dissipated. Recent developments bring to 
light the complex nature of energy exchange even in homogeneous flows 
and it seems that the concept of one-sided energy cascade has to be 
reexamined [50]. Typically, in some experiments extended areas of 
'negative' turbulent production has been revealed [51]. This is very 
pertinent to reactive flows where due to reactions, energy is absorbed 
and released in a manner that could possibly enhance the reverse cascade 
phenomena . 

The fine grained turbulence derives energy from the mean and organ- 
ized motion and dissipates some of it but returns some of the energy 
to support the generation of the larger structures. The turbulent 
kinetic energy is contained mainly in the larger eddies and only a small 
amount resides in the smaller eddies where viscous dissipation is impor- 
tant. On the other hand, for chemical reactions the species must first 
mix at the molecular level which is irore likely to occur in smaller 
eddies. Therefore spectral variation of turbulent total enthalpy 

— T 

fluctuation, H" , may be different from the variation of kinetic energy, 
— ^ 

q" . However, at present we disregard this difference for the sake of 
simplicity. 

In the determination of the dissipation rates e and we there- 
fore assume that there are two distinct ranges of wave numbers. There 
exists a range of smaller wave numbers (large eddies) in which viscous 
dissipation is negligible. In the other end of the wave number spectrum 
there is a range of large wave number (small eddies) tn which viscous 
dissipation is important. Between these two is the inertial subrange 
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which contains most of the turbulent kinetic energy q” and turbulent 
— y 

total enthalpy H" . 

In the past, investigators studied turbulent energy transfer by ap- 
proximating Reynold stresses in terms of turbulent kinetic energy, k 


au." - 

and viscous dissipation e — ) such that the characteristic time 

scale for the turbulent field was (k/e). However, it is generally recog- 
nized that turbulent interaction time scales vary in different parts of 
the energy spectrum and single scale model accounts for only the energy 
containing large eddies. Based on physical arguments and simple closure 
model Hanjalic et al. [49] developed a multiple-scale model which showed 
striking improvement in the level of agreement with experiment over that 
obtained with single scale models [49]- Essentially, they partition the 
turbulent kinetic energy spectrum E{k:) between the production region 
(k < tc-l) and transfer region (k-j < k < K 2 ) with negligible energy in the 
viscous dissipation region (< > < 2 ^ such that the total turbulent kinetic 
energy is given by Figure 40. 


k = (u^ + = 0-'' ' E(k) die E(k) d(ic) = kp + k.p 


(113) 


where ki is the wave number at which partitioning takes place. 

Energy is transported from the production (low wave number) region 

at a rate e and enters the dissipation region (high wave number) at a 
P 

rate e. The energy transfer rate in the transition region is [49]. 

In terms of the energy spectrum these dissipation terms may be written 
as 

E = 2v /"K^E(K)dk (114) 

0 

and 
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(115) 
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+ E ' Z'^*D(<) dk+ /"d(k) dk 
p So 

where D{k) would then represent the partitioned energy spectrum in terms 

of k and kT . Balance equations were modeled for the kinetic energy and 
P ' 

dissipation rates in the production and transition regions and for spec- 
tral equilibrium in very small eddies, the rate e becomes identically 
equal to . The various improvements possible with such a multiple 
scale modeling was discussed in reference [49], and specially for reac- 
tive flows this model's ability to account for different time scales of 
interaction seems promising. Beside characteristic time scales of kinet- 
ic energy transfer, in chemically reacting flows, the model must also 
allow for the characteristic time scale of chemical production rate. In 
the present formulation, this is attempted by considering a turbulent 
total enthalpy fluctuation spectrum N(k) such that 

= /" N(k) dk = N(<) dk + N(tc) dk 

0 0 1 

(11 

= (H^)p+ (H"^)s . 

This spectrum is also partitioned in a manner similar to the kinetic 
energy spectrum (Figure 40). The spectral transfer rate of enthalpy is 
then in the production region, in the transition region and 

in the dissipation region. Balance equations are then written for these 
rates along with the equations for energy dissipation rates. The basic 
modeling formulation follows Hanjalic et.al.[49] and the effect of the 
density fluctuation is modeled in a very simple manner. Since the pres- 
ent formulation is basically an integral one, the differential balance 
equations are integrated across the shear layer which results in the fol- 
lowing integral equations. 
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a) Dissipation in the Production zone: 
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where iP and jPm are the production integrals for the turbulent kinetic 
P Np 

energy and turbulent total enthalpy respectively and are given as: 
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iP-, and lP„ are the decay integrals for e and e„ respectively. 


They are given by 
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jPn and lP„ are diffusion integrals and are written as 
^i "oi 
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( 121 ) 
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Here, C„ , C„ , Cm , C„ . C and C are constants defined later. 
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jP- and are the generation terms due to mean vorticity. 
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where e , is the alternating third order tensor and 3fq = ~ 

£mk ^ 

This term vanishes for an irrotational flow and was found necessary to 
include for rotational flows [49]. 


b) Dissipation in the Transition zone: 
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here and I^m are the production integrals for the dissipation 


rates e and Em in the subrange. They are 
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and are the decay integrals and are given as 
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And finally, I^_,. and I^u are the diffusion flux integrals in 
Di Hp. 
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The constants appearing in Equations (117) to (127) are defined, 
at present, using the data available from Hanjalic et. al. [49]. They 
are given as 
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The X and are partition coefficients which define the amount 
3 f turbulent kinetic energy and turbulent total enthalpy present In the 
production zone and transition zone. Since we assumed negllgib^energy 
in the dissipation zone (< > k ^) we may thus write In terms of q" and 
2 

ull^ 
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The coefficients and Xj^ are given inputs to the problem whereas 
and are obtained as a part of the solution. 

Here X and Xn characterizes the shape of energy and total enthalpy 
q M 

spectrums respectively and and gives their degrees of spectral 
imbalance. 

The solution to the above integrals requires appropriate choices 
of shape functions. However, once the shape functions are defined, 
we can solve the complete interaction problem. Although the dissipa- 
tion rates in the dissipation zone, e and are different from the 
rates in the production and transition zones, at present, we assume 
spectral equilibrium between the transition and dissipation zones and 
take 

E = e_ 

s (132) 

^H = ^H^ 

With the assumption we can then determine the eddy diffivities € 

and €„ defined in equations (104) and (1121- 
H 

The solution of equations (117), (118), (123) and (124) along with 
the mean motion equations and turbulent fluctuation equations gives the 


42 



complete interaction problem. The problem will be completely defined 
once we can specify the shape functions for the unknown variables. Here 
we are guided by available experimental data wherever possible. However, 
there is a lack of well defined experiments specifically for reactive 
flows and so the choice of shape functions here have to be modified 
depending upon the results of numerical experiments. 

2.6.6 Shape Functions for Fine-Scale Turbulence 

The governing integral equations (67), (76), (87) and (106) for the 
turbulent correlations are solved by assuming appropriate shape func- 
tions for the correlations. Though it would be more accurate to use 
strip-wise shape functions, i.e., to do a strip-wise integration, at 
present, to keep the problem tractable and consistent with the mean 
flow formulation a single shape function is assumed for the whole shear 
layer. Furthermore, since we are assuming that all fluctuation corre- 
lations vanish as n we take the following shape functions 
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£ — = A a (x) e"’^ ; A = 1.0 

- p p 
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H"^ = N (x) e"^ 


and 


q" = Aq E (x) e ^ 


Thus o(x), C|^(x), N(x) and E(x) are the amplitude functions which are the 

basic unknowns of the problem. o(x) is the centerline magnitude of 

— 5 " 

the density fluctuation p" and similarly Cj^(x) are the centerline 

^ 

magnitude of the species mass fraction Yj^" , N(x) is the total turbulent 
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enthalpy density across the shear layer normalized by 
N(x) - J dn 
so that = }//n 

Similarly, E(x) is defined as the turbulent kinetic energy per unit length 
across the shear layer and is normalized as follows 


E(x) 




dn 


(135) 


so that Aq = 2//ir. 

We use experimental data to approximate velocity fluctuation correlations 
[52] such that 


and 
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a-j = 0.3 



(136) 


The shape functions for the dissipation rates are given by 

e. = A D.(x) e’^ 

1 * 

where i = 1, 4 denotes p, s. Hp and respectively. Again we use shape 
functions similar to the ones used for other unknown quantities. If we 
identify as the dissipation rate per unit length across the shear 
layer we can normalize equation (137) to get A^_ = l/v'iT 

The interaction problem then is basically completed by the definition 
of the various shape functions and its use in the governing integral 
equations. The reactive flow can then be solved by numerically integrating 
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the mean motion equations and the turbulent correlation equations simul- 
taneously along with the eigenvalue solution. The various constants that 
appear in the equations have to be adjusted to get the best fit with exper- 
imental data. 

The above closure method requires determination of the cross cor- 
relations of turbulent fluctuation. From existing data, for example, 

[6, 53] we see that one form of correlating data is to use correlation 
coefficients. This if and are two turbulent fluctuations we de- 
fine correlation coefficient r^^ as 


♦1' ~7~T,l/2 




(138) 


Experiments indicate that r^^ Is well defined in a shear layer and some 
guidance is available from available data [6,46]. Furthermore, except 
for the near wake area, where the interaction is highly nonlinear, we may 
consider nearly self similar profiles for r^^ across the shear layer. 
Therefore, to the first approximation we consider r^^ * specified 

as a known input to the problem. With the specification of r^^ we then 
can determine all cross correlations of the form terms of known 

variables and Numerical experiments with reacting flows [6] 

indicate that as the flow proceeds far downstream the mean flow ap- 
proaches self similarity and the correlation coefficients, r^^ tend 
to either +1 or -1. For self similar conditions, therefore, r^^ can 
be taken to be a given input to the problem. For turbulent fluctuation 
correlations of 0(3) we follow gradient formulation similar to Varma, 
et. al. [44] to achieve closure. All higher order terms are at 
present neglected for simplicity. 
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2.7 Reaction Model 


The reaction model used in this analysis is the quasi-global two-re- 
action, four-specie model due to Rogers and Chinitz [54]. This model is 
simple enough to handle numerically and has been shown to have reasonable 
agreement with experimental data. At present, this model is sufficient 
to show the effect , of finite-rate kinetics on multi -component flow. The 
reactions are given by 
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(139) 
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where k, and k. are the forward and backward rate constants and are 

f, b, 

given in reference [54] and will not be repeated here. 

The instantaneous volume rate of production of each specie is then 
determined from the reaction model (139) and are given in their non- 
dimensional form as 
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r® = E r® - — F r° 

H 2 O 2 2 '^H20 


(140) 


r° = - 2 G r° - Hr® 

OH O 2 H 2 O 
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where 
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(141) 


(142) 


On time averaging eguation (140) we get the mean rate of production 
of the i-th species, r.° appearing in equation (17). Following the fil- 
tering procedure we can separate out the contribution r^®' and r^®“ cor- 
responding to the organized structure motion and fine scale turbulent 
motion respectively. 

To the first approximation, any fluctuation in and due to 
the coherent structure and turbulence fluctuations can be neglected. 


However, when turbulence is included, the rate constants based on mean 
temperature k(T) has been shown to be significantly in error when com- 
pared to the actual mean rate constant k(T) [6]. The time averaged rate 
constants and IT. can be determined by using a probability density 

Ti 

distribution approach first suggested by Chinitz [55]. 
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In essence, a probability distribution function P(6j x, y) is con- 
sidered at any given point (x,y) to define the temperature between e and 

e + de so that 


where 6 = 


kr (6) = f k, (e) p(e) de 

^i 0 ^i 

; 0 < e < 1 


T - 

max ” "^min 


(143) 


(144) 


To solve (143) P(e) must be prescribed a priori. Chinitz [55] used a 
clipped Gaussian distribution for P(e). In free shear flows, such a dis- 
tribution seems reasonable [56] and therefore was used here also. Here, 

T is the maximum temperature attainable in the reaction model and 
max 

is the temperature in the shear layer for no reactions. A clipped 
Gaussian distribution for P(e) can be written as 

.2 

2a 


2 

P(e) = 5(0) +11^2 


where (H(l) - H(0)) is the Heavyside function given as unity in the range 
0 < 6 < 1 9 nd zero elsewhere. Furthennoret 




(146) 


*2 ■ 7 "•f' 




Using equation (145) in the definition of 0 and 6 we get [56] 

2 1/2 _ 2 -n - 

d = texp(^) - exp( '■ JJM] 


(147) 
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and 


2 1/2 2 

e'^ + 6"^ = 1^2 " ’^'3) (fifj') GxpC""^) 

2a 

2 1/2 /, ^2 „ 

-d + p)(y exp[-i^^^] - 

2a 


(148) 


6 and (= e'^ + e"^) are related to T and (T'^ + T"^) and equations 

(147) and (148) can be inverted to obtain y and o which can be used in 

equation (145) and hence finally determine (e) by equation (143). 

'i 

_ 2 
Solution of equation (24) gives T(x,n) in the shear layer. T' 

? 2 2 
(= 2lAl |T] ) can be obtained from the eigenvalue problem (55) and T" 

— 7 

can be related to H" such that 


\ [H^- u^? - i: (m.DJ+e^ )^ 

M k f (149) 

dY - ’ ^ 

■ KV " «M IF " " « 


where gradient closure is employed for turbulent correlations u"Y|^" and 
and correlations T"Y|^" and higher order correlations are neglected 
at present for simplicity. Since the solution of the mean flow equations 
is coupled to the determination of the mean reaction rate constants, 
iteration is required for the complete solution. However, for numerical 
simplicity, the mean rate constants can be determined at any given point 
(x^) in terms of the values of T and (T'^ + T"^) from the previous stream- 
wise location (><i_i) and the iteration procedure can be avoided at the 
expense of some numerical accuracy. 
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III. METHOD OF SOLUTION 


3.1. Reactive Eigenvalue Problem 

The stability equation in terms of pressure perturbation p(n), 
equation (55 ) must be solved subject to appropriate boundary condi- 
tions. Since detailed stability analysis of reactive wakes was not 
available, extensive numerical calculations were carried out to de- 
termine the eigenvalues of the problem. In essence, for a given M^ 
and frequency B and with assumed mean flow profiles, equation (49 ), 
the eigenvalue problem can be solved using Runge-Kutta Integration 
technique. The boundary conditions used are 
(i) On the wake axis (n = 0) 

For Anti symmetrical Oscillations 
P(0) = 0 

For Symmetrical Oscillations 
v(0) = p'(0) = 0 

(ii) Far from wake axis (n •+■ “) 

P(n) -► 0 

The boundary conditions are homogeneous and therefore are not 
sufficient to establish any solution of the governing stability equa 
tion (55 ) other than the trivial solution of zero. The stability 
problem must therefore be formulated as an eigenvalue problem, i.e., 
non-zero solutions which satisfy the boundary conditions exist only 
for certain combinations of the complex wave number a and complex 
phase velocity, c. Since the frequency p = oc is constraint to be 
real, then for fixed g, we need only determine and oj (and hence 
Cj^ and Cj) to obtain the required eigensolutions. 

It can be seen that an important analytical feature of 


(150) 

(151) 

(152) 
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equation (55 ) is the existence of a singularity at the point 
u = c (critical point). Though u is real, the singularity lies on 
the real axis in the complex-plane only for a neutral disturbance. 
However, for amplified and damped disturbances (Cj < 0), the sin- 
gularity lies in the complex n-plane and numerical integration a- 
long the real axis is not possible [60]. 

Essentially, integration begins at large value of n (==ng) and 
proceeds along an indented coutour around the critical point in the 
complex r)“Plane. The eigenvalues (aj^» “j) are obtained by a linear 
search procedure to satisfy the boundary conditions. The method of 
solution is basically the approach of Mack [60, 61]. The equations 
are separated into real and imaginary parts and solved simultaneously. 
Once the eigenvalues are determined all the eigenfunctions u, v, 

3, T and Y|^ can be calculated using the linear equations (see Appen- 
dix A). For neutral eigenvalues, i.e.,c^ = 0,a linear search is car- 
ried out in Cj^ and because aj = 0 for the spatially developing 
disturbance studied here. 

3.2. Interaction Between Coherent Structure and Mean Motion in 

Reactive Flow 

When all turbulent correlations are neglected in the govern- 
ing mean flow equations, the integral equations show the balance 
between the mean motion and coherent structure motion. The inte- 
grals appearing in the equations can be calculated in terms of the 
eigenfunctions (Section 3.1 ) by Simpson's Rule. To avoid inte- 
gration to large values of n, n^ is taken to be = 5.0, since for 
n>5 the mean flow approaches free stream and the integrals in the 

range n < n ^ “ can be shown to be a function of the value of 
^ 6 
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eigenfunction at n = n^. The governing equations for the mean flow 

2 

variables V_, T and and coherent structure amplitude |A] can be 

C C 

written in a general functional form as 


f.(4>^ ;Re , i = l,7 (153) 

where denotes the unknown variables of the flow, i.e., V^. T^.; Y|^ 

(k=l,...4) and [A] . These equations can then be solved in the stream- 
wise direction from some initial conditions at, say, x = x^ such that 
V,(x„) = Vc ,T^(Xo)-T and |A|2 (x„)=|A|/. Since the 

0 0 C Cq 

eigenvalue problem is a function of <f>^, at each streamwise location the 
eigenfunctions and the integrals are evaluated to proceed with the solu- 
tion of equation (153). The various integrals appearing in the equa- 
tions are given in Appendix B. 

3.3. Fine-Scale Turbulence in Reacting Flow 

The integral equations for turbulent-correlations, equations ( 67 ) , 

( 76 ), ( 87 ) , and (106), can also be written in a similar form to equa- 
tion (153) and solved subject to proper choice of initial conditions for 
the shape functions (133). In all, there are 11 equations that must be 
solved along with the mean motion equations described in Section 2.3. How 
ever, at present the turbulence model has not been included in the numeri- 
cal calculation of the interaction problem because of lack of comparison 
with experimental data, specially for reactive wake. Instead, a numerical 
parametric study of the turbulent integral equations was carried out to 
determine the effect of the variation of the constants and correlation 
coefficients R^^ appearing in the equations. For simplicitly, the mean 
flow was assumed to be self similar in transformed coordinate n and in- 
teraction between the coherent structure and turbulence was neglected. 
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For self similar conditions, all cross correlation coefficients 
tends to either -1 or +1 and therefore from past studies [6], was 
given a constant value across the shear layer. Effect of the variation 
of R on the development of the turbulence field is studied parametri- 
cally. Effect of adjusting the various constants (i.e.,Ep , Sj^, S^) 

was not studied in detail due to lack of experimental data for the physi- 
cal reaction model employed here. Numerical calculations were carried 
out using the Runge-Kutta integration technique. 
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IV. RESULTS AND DISCUSSION 


In this chapter the results of the numerical calculations are pre- 
sented and the various facets of the interaction problem are discussed. 

Since the basic thrust of this research has been to develop a model to 
study the interaction stability of a coherent structure in a supersonic 
reactive flow, detailed analysis of the stability of a disturbance in 
reactive flow has been carried out. There is a considerable lack of ex- 
perimental data on the reactive wake phenomena, specially in the H 2 - O 2 
combustion system studied here. Therefore, no direct comparison was pos- 
sible for the development of a coherent structure in reactive free shear 
layer. However, experimental data [27,28,40] is available for nonreactive 
hypersonic wake of a flat plate and a comparison is possible. The reac- 
tion model is at present not of great importance to study the growth and 
decay of coherent structure except so far as to include finite-rate kine- 
tics in the interaction problem. The presence of finite-rate chemistry 
results in large species production rates for the reaction model used 
[54] and therefore imposes stringent stability restrictions on the expli- 
cit numerical calculations carried out here. The step size required for 
stability of the numerical calculations of the reactive interaction pro- 
blem is so small that detailed computations for reactive cases was compu- 
tationally prohibitive and therefore only a representative set of solutions 
are presented. It may be possible to solve this problem by an implicit 
numerical technique which could result in detailed solutions within rea- 
sonable computational time. At present such a numerical solution is be- 
yond the scope of this research and must await future study. 

The turbulence model developed in this study is at present untested, 
and therefore its inclusion into the interaction between coherent structure 
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and mean motion cannot be justified. Instead, the turbulence model is 
parametrically studied for a simplified mean flow and with the assumption 
that there is no interaction between coherent structure motion and fine- 
scale turbulence. This enables the fine-scale turbulent motion to be- 
come uncoupled from the development of the instability wave and with pro- 
per choice for the mean flow profiles, the turbulent integral equations can 
be numerically solved subject to given initial conditions. The various 
constants and correlation coefficients appearing in the turbulent 
equations have been parametrically varied to study their effect on the 
spatial development of fine-scale turbulence in reactive wakes. De- 
tailed comparison with experimental data for a free shear layers under- 
going ^^2 ~ ^2 again not possible due to the lack of experimen- 

tal data on the turbulence field of this combustion system. The success 
of any turbulence model capable of handling finite-rate kinetic effects 
will depend upon the availability of reasonable comparison between nu- 
merical predictions and experimental data and since this is not possible 
at present, the turbulence model developed here Is studied only as an in- 
dependent problem and was not included in the interaction analysis. 

The results presented here therefore falls in three distinct parts. 

In the first section the results of the detailed study of the stability 
of reactive wakes is presented. In the second section, the interaction 
between the coherent structure motion and the mean flow motion is numeri- 
cally studied. And, finally, the numerical parametric study of the tur- 
bulence model is presented for a self similar mean flow field. Effect 

of varying the cross correlation coefficients is studied for both reactive 
and non reactive cases. 
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4.1. Stability of Reactive Laminar Wakes 

The eigenvalue problem developed for reactive flow, equation ( 55 ) , 
is solved subject to the boundary conditions (150), (151) and (152) and 
the eigenvalues a(ap, oj) are numerically calculated. With proper choice 
of the the mean flow variables (i.e., V^, and Yj^ ) and for a given 

free stream Mach number M , the stability equations can be integrated 

00 

as described in Section 3.1. Solution proceeds in the complex n-plane 
(Figure 2) and therefore the equations and the boundary conditions must 
be separated in their real and imaginary parts and solved simultaneously. 
With a reasonably good choice of initial guess for and aj for a given 
real frequency 3, the nonreactive case converged in 4-6 iterations. How- 
ever, the corresponding reactive cases took 6-9 iterations and the compu- 
tational time also increased proportionately. A convergence criteria of 
10"® was used for both cases. 

It was seen from the computations that the eigenvalues for the re- 
active case depended significantly on the finite-rate kinetics used in 
this model and therefore it is expected that different rate kinetics me- 
chanism would result in different eigensolutions. However, it was seen 
that the general trend of solution remains the same and so some generalized 
conclusions can be drawn from the study of this reactive stability pro- 
blem. Stability calculations such as this gives some insight into the 
phenomenon of transition to turbulence in reactive flows which is an 

area of research of great interest. 

Lees and Lin [43] and Lees and Gold [57] have shown for nonreactive 
flows, that neutral subsonic disturbances exist only when a generalized 
mean density-vorticity product has an extremum, and then only when the 
phase velocity corresponding to this location is subsonic relative to an 
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observer fixed in fluid. This condition is also sufficient (but not neces 
sary) for the existence of amplified subsonic disturbances. It has been 
shown that the relative velocity 1 - u(x,0) governs the stability of the 
wake and so the proper coordinates for the stability problem is a system 
fixed "in the fluid" even though the mean flow is calculated in a coor- 
dinate system fixed in the body. A disturbance is classified "subsonic," 
"sonic" or "supersonic" according to whether the relative propagation 
velocity is less than, equal to, or greater than the ambient sound speed. 

The local relative Mach number of the wave front is given by [57] 
M,(w-c,) 


M = 


(154) 


/T 


where M = M V is the relative Mach number. A subsonic disturbance must 
r c 


therefore satisfy the condition 


(155) 


It was also shown for nonreactive case that a neutral subsonic distur- 
bance exists only for certain mean velocity temperature profiles and -Cp 
is uniquely determined by these profiles. The stability restriction re- 
sults in the criteria that -Cj^ must lie in the interval [0,1] for a 
self-excited disturbance to exist. For subsonic disturbance it results 
in the requirement that [57,58] 

0 < -Cp < i;^ (156) 

It has been possible to establish some criteria for the existence 
of neutral subsonic disturbances [57] in compressible nonreactive wakes 
but it seems impossible to extend such requirements for the reactive 
case. For two-dimensional compressible nonreactive wakes, the condition 
that uniquely determines the existence of neutral disturbance is given 


by 
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(^) 


and 


T c, 


"r ' "R. 


(157) 


where prime denotes differentiation with respect to n* Using the mean 
flow profiles and equation (157) results in 




-n. 


- ■> 

(-^) 


and 




-n. 


(158) 


(159) 


where is the critical point where w = Cj^. 

This indicates that for nonreactive case, the neutral solution 

phase velocity Cj^ is a function of only and is independent of the 
c 

relative Mach number. For these values of c„ , the neutral inviscid ei- 

c 

genvalues for anti -symmetric disturbance (a^) can be determined as a 

function of T^and (Figures 3,4) . The solution shown in Figure 4 

clearly indicates that whereas for no reactions in the flow, the criteria 

equation (157) is always satisfied, for finite-rate kinetics however, 

the neutral phase velocity, -c. , is no longer a function of T^. only but 

c 

instead depends on both and M^. This interesting deviation from the 
inviscid stability criteria, equation (157), indicates that the require- 
ment given by equation (157) is necessary but not sufficient to theore- 
tically define stability bounds for reactive wakes. At present, an ex- 
tension of the nonreactive stability criteria to handle reactive stability 


seem very difficult and has not been attempted. 

In Figure 3 results of the neutral eigenvalue solutions for nonre- 
active flow is presented. For increasing temperature excess T^., -Cp^ in- 
creases towards zero, i.e., the phase velocity approaches free-stream 


( 
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velocity. The critical point where w= Cj^ increases with T^.. Excellent 

c 

agreement with the Lees and Gold [573 solution also obtained (Figures). 

The eigenvalue decreases and almost becomes constant with increase in 

T for the case of m 3 = 0. The dependence of the phase velocity on T^ 
c r 

and is shown in Figure 4 where it can be seen that the reactive case 

deviates significantly from the nonreactive case. Clearly, the numerical 

solution shows that a minimum value for -c„ for each value of T^ is ob- 

c 

tained as a function of {* neutral eigenvalues calcu- 

lations a value of V^. = 0.4 was maintained at all times. It is also noted 

from Figure 4 that the -Cn bottoms out much more sharply for lower T^. 

^c 

The actual phase velocity, in the coordinates fixed with the body is 
given by Cj^=l + V^.Cj^ and therefore depends implicitly on the velocity de- 
fect in the wake axis. In Figure 5 the solutions for the neutral eigen- 
values a as a function of M_ and T_ is presented. For fixed value of 
C r c 

T the wave number decreases with increasing M as predicted by Lin 
c c ‘ 

[593 and computed by Lees and Gold [573. Again excellent agreement with 
the solutions of reference [573 is obtained for the nonreactive case. 

For reactive case, varies with in a similar manner as the nonreac- 
tive case. However, the wave numbers are much larger and the range of 
relative Mach number over which the neutral (and adjacent amplified) dis- 
turbances can exist also increases. This kind of increase is also noted 

with the increase in temperature excess T^.. 

Since we are interested in spatially developing waves, a comparison 
of the variation of frequency e with relative Mach number M^ is presented 
in Figure 6 . For fixed T^., e decreases with M^ for both reactive and 
nonreactive cases. The frequency at which neutral disturbance can exist 
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increases with increase in T^. and is also greater for the reactive case. 
The leveling of frequency B for T^=4.0 during reactions is due to the 

variation of -Cj^ observed in Figure 4. 

The boundary conditions, equations (150)- (152) allow for both sym- 
metrical and anti symmetrical disturbances in the wake and therefore nu- 
merical solutions were obtained for both the disturbances. It has been 
determined that anti synmetri cal disturbances are more unstable in laminar 
wake transition [29] and similar trend is also observed during the pre- 
sent computation. Figures ( 7 ) to (10) give the eigenvalues computed 
for symmetrical disturbances and Figures (11) to (14) give the corres- 
ponding anti symmetrical disturbance solutions. The complex phase velo- 
city C(Cf^,Cj) shown in these figures is the velocity with reference to 
the coordinate system fixed to the body and is defined as C=l+V^Rec 

where c (Cn, c.) is the relative phase velocity defined in the local 
K 1 

coordinate system (equation (55)). 

In general, the reference conditions used for the present compu- 
tations unless otherwise specified are 

T = 1500°K , Y,, = 0.05, 0.90, 0.025, 0.025 

L = 0.1m M = 2.0 and B = 0.2 

CD » 

and the centerline values for the mean profile were taken to be 

V =0.4 , T =4.0 , Yi, =(-0.25 ,-0.25,3.0 , 6.5) 

^ c 

In Figures ( 7 ) to ( 9 ) the reference values for species Y,^ was 

taken to be (0.3 , 0.65 , 0.025, 0.025) which is a fuel rich mixture for 
this present reaction model. 

In Figure 7 we see that for nonreactive flow, with b < 0.24 the 
solutions are amplified (Cj > 0) and with B > 0.24 the solutions are 
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damped (Cj < 0). At e = 0.24 we have a neutral eigenvalue solution. 
However, when the reaction is turned on, all the solutions are amplified 
and also the degree of amplification is much higher. The maximum ampli- 
fication occurs around g = 0.092 for nonreactive case but is around 
g = 0.081 for the corresponding reactive case. This seems to indicate 
that the disturbance is more likely to become unstable for the reactive 
case. The phase velocity monotomically increases with g for both 

cases. However, is lower when chemical reactions are going on. This 
implies that for reactive flow, though the disturbance gets highly ampli- 
fied it is travelling slower relative to freestream. This is consistent 
with experimental observation by Ganji and Sawyer [21] of organized 
structure motion in a propane-air combustion flow. 

In Figure 8 we show the variation of the eigenvalues as a function 

of freestream Mach number M . A neutral eigenvalue is observed near 

00 

M =2.375 for nonreactive case which disappears when the reaction is 

oo 

turned on. We notice that for M <2.375 the solutions are amplified 

00 

(C, >0) and for M > 2.375 the solutions are damped. This is very inter- 
esting because in a very early paper concerning stability of jet flows, 
Pai [62] had theoretically predicted that there is a critical Mach num- 
ber above which jet type flow will be stable with respect to all small 
symmetrical disturbances in the inviscid fluid. He further showed that 
for the case of no energy transfer the critical Mach number is about 2.5. 
Our present numerical results for symmetrical disturbances shows that 
the flow is stable for M^> 2.375. This close agreement between earlier 
theory [62] and present analysis indicates the validity of the present 
numerical results. Another interesting observation in Figure 8 is that 
when the reaction is turned on there is no neutral eigensolution and the 
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results are completely amplified. Moreover, the reactive solutions 
seems to be weakly dependent on A possible explanation of this re- 
sult is that since the release of chemical energy due to reactions is 
independent of Mach number, the variation of M^ will probably have limi- 
ted effect on the eigenvalues. The phase velocity is still smaller 
when reactions is turned on. 

Figure 9 shows the variation of the eigenvalues as a function of 
velocity defect V^. There is significant differences between the non- 
reactive and reactive cases. For nonreactive flow the rate of amplifi- 
cation is maximum around V^= 0.505. However for reactive cases calcula- 
ted here, the disturbances amplification increases continuously with in- 
crease in V^. The phase velocity is again smaller for the reactive 
case. The results indicate that very close to the rear edge of the flat 
plate where is large, chemical reactions make the disturbance grow 
rapidly but on proceeding downstream into the wake where the interaction 
weakens and decreases, the amplification rate decreases. The phase 
velocity however, increases as is decreased and approaches free- 
stream velocity far downstream into the wake. 

We further note that for V^<0.45, the amplification rate is lower 
for reactive case as compared to reactive case. This indicates that the 
instability wave is not amplified for all cases considered here and that 
under some initial conditions reactive wave will be less amplified when 
compared with the nonreactive case. It was determined during the pre- 
sent calculations that the wave amplification is quite dependent on 
shape function assumed. The phase velocity, however, is always lower 
for reactive cases. In Figure 10 we plot the variation of the eigen- 
values as a function of temperature increment T^. Here for small value 
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of we see that Cj is larger for nonreactive case as compared with 
the reactive case. The phase velocity is still lower for reactive 
case. There is a maximum amplification occuring around T^=2.4 for 
nonreactive case but is no longer present for the reactive case. This 
figure shows that with the low value of considered here, and for reac- 
tive case with 5 %H 2 ,the wave is amplified less as compared with nonreactive 
case but the amplification continuously increases with increase in T^. 

With higher value of V^, the amplification again becomes higher for reac- 
tive case as compared with nonreactive case. The amplification rates 
show similar trends as observed experimentally by Ganji and Sawyers 
[ 21 ]. 

In Figure 11 we present representative calculations of the eigen- 
values for the case of antisymmetric boundary conditions, equation (151), 
as a function of 3. In previous wake studies it was determined that 
antisymmetric boundary conditions resulted in more unstable wave motion 
and it has also been experimentally confirmed that in wake flows the 
wave structure correspond more closely to the antisymmetric solution [27 ,28 
29]. As can be seen in the figure, the amplification (Cj>0) is much larger 
for reactive case as compared to nonreactive case. Furthermore, when 
compared with the symmetric solutions in Figure 7 we see that amplifi- 
cation is greater for the antisymmetric solution. All the eigenvalues 
are amplified solutions and no neutral solution exist unlike the case 
in Figure 7 . For nonreactive case there is a maximum amplification 
occuring around 3 = 0.17 which is greater than 0.092 for the corres- 
ponding symmetrical case. The wave speed Cp is still smaller for reac- 
tive case as compared with the nonreactive case. 

Figure 12 presents the eigenvalues as a function of for 


63 



anti synmetrl cal disturbances. The amplification continuously decreases 
with increase in M for no reactions and shows a maximum for very low 
Mach number. For reactions the amplification shows a maximum around 
M^=3.0 and is generally much larger than the nonreactive solution. 
Compared with the corresponding symmetrical disturbance solutions 
(Figure 8), anti symmetrical disturbance is highly dependent on the 
freestream Mach number. Figure 8 also shows the stabilizing effect of 
increasing M^ as observed experimentally [27,40]. 

Figure 13 shows the eigenvalues as a function of velocity defect 
V^. In general the amplification rate is higher for reactive case and 
the overall behavior is quite similar to the corresponding symmetrical 
case (Figure 9 ). However, no maximum is directly observed in ampli- 
fication case within the range calculated except perhaps the tendency 
to level off near V^ = 0.8. Figure 14 gives the variation of the eigen- 
values as a function of T . The trend shown here is quite different 
from the symmetric case. The amplification for nonreactive case de- 
creases with increase in temperature excess as observed by Lees and 
Gold [57] and the wave number also decreases with T^. For reactive 
case however, the amplification rate shows a distinct maximum around 
Tj.= 1.8 and the wave number variation also shows a maximum. This in- 
dicates that the wavelength (= l/oj^) first decreases and then increases 
with increasing temperature excess. The phase velocity Cj^ is consis- 
tently lower for the reactive case as determined by Ganji and Sawyer 
[ 21 ]. 

Finally, Figures 15 and 16 show the characteristic variation of 

A 

the pressure eigenvalue amplitude |p| and phase ifip across the shear 
layer for symmetric and antisymmetric disturbance respectively. In 
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general, the amplitude is much larger for the reactive case in both 
figures. The phase changes across the shear layer shows significant 
difference between reactive and nonreactive flows. 

Finally in conclusion the study of the stability of reactive 
wakes has shown remarkable differences with the classical nonreactive 
stability analysis. Though small disturbance analysis cannot deal with 
transition, it however has been shown [39] that linear stability theory 
is a useful guide to the understanding of the main parameters and phe- 
nomena governing transition. The growth and decay of amplification 
rates along the wake axis can be used to give a good indication of non- 
linear effects present in the wake. The excellent agreement of lami- 
nar stability theory with previous studies for nonreactive cases [57] 
indicates the applicability of this theory to study the development of 
disturbances in the flow. The present calculation of the eigenvalues 
shows good agreement with past theoretical and experimental analysis. 
For reactive flows, it has been shown that the neutral phase velocity 
is no longer independent of relative Mach number M^. The condition 
for neutral stability equation (157), is no longer sufficient for 
reactive flows. Furthermore for the amplified solutions calculated 
here, qualitative agreement has been obtained with the experimental ob- 
servations of Ganji and Sawyer [21]. In essence, the reacting wakes 
are more unstable than nonreactive wakes. The amplification of the 
wave is in general higher for reactive wakes but the wave propagation 
velocity is consistently lower. The antisymmetric disturbance is more 
unstable than the symmetrical disturbance in general for both reactive 
and nonreactive cases which again agrees with experimentally observed 
data in wake flows. The calculations further show that the choice of 
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the mean flow velocity, temperature and mass fraction profiles governs 
to a great extent the spatial stability of laminar reactive and non reac- 
tive wakes. 

Since the reactive eigenvalue problem formulated in Appendix A is a 
very general one it was decided to do some representative calculations of 
reacting flow fields other than a wake. Essentially, the mean flow as- 
sumption has to be modified to account for different physical geometry. 

For example, a heated symmetric jet has a mean velocity flow profile that 

_ 2 

can be given in terms of similarity variable n as u = sech n [67]. For 
heated jet, the temperature profile can be still taken as near Gaussian, 
equation (40). For symmetrical disturbances, the stability problem was 
solved for a characteristic frequency 6 = 0.1 and presented in Figures 17 
and 18. Figure 17 shows the variation of the amplitudes of pressure 
eigenfunction [p| and vertical velocity eigenfunction |v[ across the jet 
for reactive and non reactive conditions. Comparing with the symmetric 
solutions for the wake. Figure 15, v;e see remarkable differences that do 
not necessarily depend upon the differences in reference conditions alone 
The vertical velocity amplitude |v| shows a maximum for non reactive case 
around n 5. However, for reactive cases no such maximum is observed 
and the amplitude is lower than the non reactive case. The eigenvalues 
indicate that amplification is much higher for reactive case and the 
wavelength (~ is lower. The corresponding phase changes of the 

pressure and vertical velocity eigenfunction across the jet is shown in 
Figure 18. The pressure phase shows a single phase change at large 
value of n for both reactive and non reactive cases. On comparison with 
the corresponding wake solution. Figure 15, it can be seen that reactive 
wake has two phase changes. The phase changes for the vertical velocity 
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eigenfunction 41 ^ is also shown in Figure 18. In general, the trends in 
phase changes are similar for both non reactive and reactive jets. De- 
tailed calculations of reactive jet stability has not been carried out 
at present. However, it has been shown that the generalized stability 
problem formulated here is quite versatile and is capable of handling 
different physical problems. Other reactive models can also be handled 
by proper modifications to the finite-rate kinetic terms v* and F in 
equation (55). Reactive flow stability of other shear flows like boun- 
dary layers and jets would also be of great interest due to the exis- 
tence of large scale structure in these shear flows [25 ]. Such stabili- 
ty calculations have, however, not been carried out at present and must 
await future study. 

The study of wave development and the onset of nonlinearity leading 
to transition has been extensively studied and the general results show 
both qualitative and quantitative agreement with past experimental and 
theoretical data. The present study results in some basic insight into 
the instability wave phenomena in reacting shear layers. By modeling 
coherent structure as an instability wave the results of such an analysis 
can be used to study the effect of finite-rate kinetics on the develop- 
ment of organized motion in free shear layers. However, the spatial 
development of such instability waves is coupled to the variation of the 
mean flow and therefore their spatial development has to be studied 
simultaneously. Such an analysis is the topic of study in the next 
section. 

4.2 Coherent Structure Interaction in Reactive Wakes 

The development of the interaction between the mean flow and the 
coherent structure can be studied by integrating the governing equations 
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for given values of Re, Cp and The integration proceeds down- 
stream from an initial location x=x„ with initial values V , T , Yj^ 

° ‘'o 0 Cq 

and |A1 q^. Experimental data [27-29, 40] can be utilized to obtain ini- 
tial values for the velocity defect and temperature excess for non reac- 
tive wakes. However, the initial value for the amplitude of the coherent 
structure has not been measured experimentally and so the value used must 
be such that it can simulate experimental conditions. There is also a 
considerable lack of experimental data on reactive wakes and so the ini- 
tial values used are the same as that for non reactive wakes and a com- 
parison of the interaction development is studied between the two cases. 

For non reactive flow, with uniform concentration gradients, i.e., Y|^ -0, 

k=l,4, the interaction problem is numerically calculated for the experi- 
mental data available [27, 40]. Using the initial conditions used by 
Liu and Gururaj [39] the solutions are obtained for a hypersonic wake 
[27]. The calculations correspond to the transition behind an adiabatic 
flat plate. The initial location is chosen such that it is upstream of 
the nonlinear region. Thus the interaction develops from the linear re- 
gion into the nonlinear region. In experimental studies the transition 
from laminar to turbulent arises naturally from body generated noise or 
other forms of disturbances. In the present numerical calculations the 
disturbance arises due to the initial value given for the amplitude 

and its growth and decay downstream occurs due to the interaction 
between the mean flow and the coherent structure. 

Figures 19 and 20 show, respectively, the variation of the center- 
line temperature excess and centerline Mach number along the wake axis. 

The present calculations show excellent agreement with the experimental 
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data of Batt and Kubota [27] and the theoretical work of Liu and Gururaj 
[39]. In the initial growth region the disturbance amplitude is small 
and both the temperature defect and Mach number M^ follow essentially 
the linear growth and is identical to the 'laminar' case, \A\^ =0. However, 
as the interaction proceeds downstream the amplitude rapidly increases 
and the solution deviates from the 'laminar' solution and the temperature 
and velocity decay is more rapid. The solution follows very closely the 
experimental data [27]. The decay is more rapid for higher values of 
jAlp and it was seen that the higher value of (=2x10 ) shows a 

better comparison with the experimental data [27]. Since detailed cal- 
culations of the non reactive wake has already been investigated [29, 38, 
39, 63], such calculations will not be repeated here in detail. Only 
some representative solutions for non reactive cases will be presented 
here. 

It was found during the preliminary calculations of the reactive 

wake that the reaction model used here [54] results in large production 

rates and therefore imposes stringent stability requirements on the inte- 

gration of equation (153). In fact the complete reaction model, equation 

-9 -10 

(139), results in integration step size of the order of 10 ~ 10 that 

is computationally prohibitive. In their calculations also, Rogers and 
Chinitz [54] also observed large instantaneous production of the hydroxyl 
radical thereby requiring very small integration steps. These large pro- 
duction rates were identified with the first reaction in equation (139) 
and they found typical chemical time for OH production of the order of 

S 0 conds [54]. They further estimated that the first reaction, i.e., 
^2 '*’^2^ 20H, may be in chemical equilibrium during much of the time 
required for a typical integration step that is computationally reasonable. 
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Due to the limitations imposed by the available computing system, 
it was decided that in order to alleviate the computational problem 
caused by large production rates, a simplified procedure would be adopted 
that would allow finite-rate calculations within reasonable computational 
time. The present formulation already assumes premixed flow and since 
we are only interested in presenting a formulation that can account for 
coherent structure motion in reactive flow, to the first approximation, 
we assume that the first reaction is in equilibrium and no production of 
OH radical occurs. The second reaction which is essentially the produc- 
tion of water vapor is a much slower reaction and finite-rate calcula- 
tions based on this reaction is computationally possible. Therefore, 
all reactive solutions presented here will assume the above mentioned 
simplifications. This assumption results in no production of OH radical 
and also no consumption of ^ 2 - premixed model used here did not in- 

clude as an inert specie in contrast to the model in reference [54]. 
By assuming equilibrium for the first reaction, the presence of oxygen 
essentially assumes an inert form for the second reaction since no con- 
sumption of O 2 occurs. The results presented here, therefore, account 
for the finite-rate production of water vapor at the expense of the hy- 
drogen and hydroxyl radical in the premixed flow. Finite-rate calcula- 
tions of the second reaction still results in step sizes of the order of 
1 q"4 ~iO”^ and therefore only representative calculations are presented 
due to computational system limitations. 

It must be noted here that the reactive eigenvalue problem depends 
significantly on the variation of mean flow profiles and must be solved 
at each streamwise location which results in increased computational 
costs. For numerical efficiency it is however possible to store an 
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eigenvalue map, whereby the eigenvalues (a, c), the eigenfunctions 

and their integrals can be stored as a function of V^,, and Y|^ (k= 1,4). 

c 

At each streamwise location this map or table can then be referred to 
and the interaction integrals (Appendix B) determined by interpolation. 

At present such numerical optimizations has not been carried out due to 
storage logistics. 

It was also observed during the calculations that the mean flow 
shape functions assumed may not necessarily be applicable to reactive 
flow situations. Even though the mean velocity and temperature profiles 
have been experimentally observed to be nearly Gaussian [65] no detailed 
measurements are available for the specie concentration profiles in super- 
sonic reactive wakes. Concentration profiles in reactive jets [22] indi- 
cate that the profiles do deviate significantly from Gaussian distribu- 
tion [23]. Since the integral formulation used here requires a priori 
assumptions of the mean flow shape functions, proper choice of the pro- 
files is necessary for accurate solutions. Here, for numerical simplicity 
the mean flow shape functions for all variables were assumed to be Gaus- 
sian as defined in equation (49). This simplification allowed the reduc- 
tion of all the integro-differential equations to the form in equation 
(153). During the course of the numerical calculations, it was deter- 
mined that this simplification resulted in reactive solutions that were 
highly dependent on the mean flow specie mass fraction shape function 
assumptions. Representative mass fraction profiles are shown in Figure 
21 (a-d) for a supersonic reactive wake. To see the effect of large 
concentration gradients near the wake centerline, the mass fraction values 

at n=0 were taken to be (= -0.25, -0.25, 3.0, 6.5). Figures 21a 

^c _ 

and b show the mean flow mass fraction profiles for hydrogen, Yy and 

"2 
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oxygen, at various axial locations. We see that hydrogen levels are 
reduced as the flow progresses downstream. Similar trend is observed in 
Figure 21c where the hydroxyl radical mass fraction profiles are presen- 
ted. The chemical reaction used here results in the depletion of H 2 and 
OH as the water vapor is formed. This is evident from 21d where 

profiles increase with axial distance. Comparing with the non reactive 

profiles also plotted in Figure 21, we see that all the non reactive 

profiles tend to uniform conditions with axial development. However, 

for reactive conditions, only Y^, profiles still approach uniformity 

^2 

which is directly due to the fact that the first reaction is assumed to 
be in chemical equilibrium and oxygen molecules act as an inert specie 
for the second reaction. These figures also indicate that OH radical is 
consumed more rapidly than the fuel (H 2 ). This is again consistent with 
the fact that OH is a very unstable radical and at high temperatures re- 
acts very fast to produce H 2 O vapor which is more stable. The numerical 
computation therefore had to take into account the large rate of deple- 
tion of OH as the restricting factor leading to small integration step 
sizes. When all OH radical mass fraction was used up at the centerline 
(i.e., Yqi^ — »■ -1) the numerical computations were terminated. This was 
necessary to avoid negative mass fractions in the flow field. This limi- 
tation of the present calculations can be related to the choice of mean 
flow specie shape functions and also to the fact that molecular diffusion 
was neglected in the original formulation. Another reason for this be- 
havior is the basic assumption of premixed conditions which has been 
shown to be more complicated than non-premixed flames [68]. Furthermore, 
on comparing the production rates in supersonic H 2~02 combustion systems 
with other hydrocarbon combustion systems [21] it was seen that the 
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hydrogen-oxygen combustion rates are much higher leading to the point of 
view that the combustion system used here is 'kinetically' dominated. 

However, all the above mentioned restrictions in this formulation 
can be removed by solving with more appropriate shape functions for the 
specie and adding diffusion terms, etc. Of course, finite difference 
calculations would be more accurate instead of integral considerations 
for reactive flow fields, but step size restrictions would be more 
stringent in such calculations. 

Figures 22 and 23 present respectively the centerline variation of 
the mean flow temperature excess and velocity defect for a M^=3 turbu- 
lent non reacting wake experimentally studied by Demetriades [40, 64]. 

To account for turbulence simple eddy viscosity models in Spalding [69] 
(Model I) and Schetz [70] (Model II) (see Appendix C) were utilized. 

The results indicate that the decay in temperature and velocity is only 
slightly larger than the case for no turbulence. Coherent structure 
development as discussed above cannot account for turbulence and the 
results indicate that the decay observed due to the interaction between 
the large structure and the mean flow deviates significantly from the 
turbulent data of reference [40, 64]. Better comparison is observed 
in the earlier part of the wake where the amplitude [A| is increasing. 
Further downstream, the amplitude decays and the coupling between mean 
flow and coherent structure motion weakens. This growth and decay be- 
havior of the amplitude can be seen in Figure 24 which also shows the 
balance between the various terms appearing in the coherent structure 
kinetic energy conservation. Here the averaged disturbance kinetic 
energy density C(x) is given by 

C(x) = |A1^ = 1 f ”(101^ + lavl2)dn ■ (11 
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which is due to the normalizing condition used in equation ( 59). We 
further note that the development of the amplitude is essentially 
governed by the Reynold stress production due to the interaction 
between mean flow and coherent structure. All other terms contributing 
to the conservation mechanism in equation (58) are relatively smaller 
than Very similar results were also obtained by Lin and Gururaj 

[39] and Ko et. al . [38]. Inclusion of turbulence into the governing 
equations results in lower values for the energy density C(x), stress 
production and coherent structure dissipation However, the 

mean flow dissipation I., also shown here, is larger for the turbulent 

<p 

case as expected. The mean flow dissipation continuously decreases 
with axial development since the mean flow velocity gradient (-^) de- 
creases. The coherent structure dissipation first increases with the 
increase in the amplitude |A| and reaches a maximum value at an axial 
location that closely corresponds to the amplitude maximum. Further 
downstream as the amplitude decays the dissipation, also decreases 
as v/as observed in other calculations [39]. 

The calculations for non reactive cases show only limited effects 
of the inclusion of the eddy viscosity model. The more complete turbu- 
lence model formulated in Chapter II was not used due to lack of time 
to do extensive numerical experiments. Since the present analysis is 
directed towards an investigation of the interaction between mean flow 
and coherent structure in turbulent reactive flows, the eddy viscosity 
model was included to account for some turbulent effects. Due to the 
lack of experimental data on supersonic turbulent reactive wakes, no 
detailed analysis was undertaken to investigate the effect of various 
turbulence models. However, all further results presented here for 


74 



reactive cases were calculated with inclusion of the algebraic model for 
eddy viscosity (Model II). 

4.2.1. Effect of Initial Amplitude and Oscillation Frequency 

Figures 25-28 present the results for the reactive wake for various 

2 

initial values of the amplitude, jA|^ . Figure 25 gives the centerline 

temperature increase as a function of axial distance. Since the second 

reaction is exothermic, temperature is increased as water vapor is formed. 

2 

We see that as the initial amplitude |A|q is increased, the temperature 
also increases. In contrast, non reactive flow temperature excess de- 
creases with increase in initial amplitude. From the calculation it 
became evident that the increase in amplitude results in increased con- 
tribution from the coherent structure development to the mean flow ther- 
mal energy conservation. Essentially, the streamwise gradient of total 

c 

enthalpy diffusion flux integral, - and dissipation I^ are always 
positive and increases with increasing initial amplitude which directly 
contributes to the growth of temperature excess. In constrast, the decay 
of mean flow velocity defect (Figure 26) increases with increase in ini- 
tial amplitude. This is due to the fact that the Reynold stress produc- 
tion due to coherent structure, increases with increase in initial 

amplitude and therefore increases the velocity defect decay (equation 
(36 )). The mean flow dissipation, I^, (equation (46 )) is not affected 
by changes in initial amplitude. Figures 27a and 27b give the evolution 
of the mean flow centerline specie mass fractions under the same condi- 
tions. Increase in initial amplitude increases the decay of hydrogen 
and hydroxyl centerline values (Figure 27a) and also increases the in- 
crease of water vapor (Figure 27b). However, it is clear that in this 
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kinetically controlled combustion studied here, the effect of changing 
the initial amplitude is only marginal. Essentially, this is due to 
the fact that the combustion results in large production rates which are 
almost independent of the coherent structure development. Furthermore, 
tliough the growth of the amplitude increases with increasing initial 
values (Figure 28), the maximum is significantly lower when compared to 
the non reactive case. The maximum occurs around x/L —13.2 as compared 
to x/L 30 for the non reactive case (Figure 24). One of the interes- 
ting results of such calculations is that the growth rate of coherent 
structure is significantly lower in reactive flow and also the maximum 
attained is lower. Evidently, this reaction mechanism inhibits the 
growth of coherent structure, a result consistent with the observation 
of Ganji and Sawyer [21]. Though the growth and decay behavior of the 
amplitude is similar in both non reactive and reactive flows, their 
magnitudes are quite different. This essentially seems to indicate the 
growth rate of coherent structure amplitude depends quite significantly 
on the reaction mechanism. In fact, the present calculations show that 
due to the large production rates associated with the reaction, the 
coupling between mean flow and coherent structure motion is considerably 
weakened and the mean flow field variation, especially the thermal and 
concentration fields are not significantly affected by the coherent mo- 
tion. The coherent structure motion however is modified due to the pre- 
sence of finite-rate kinetics. This behavior can again be identified 
with the dominance of the finite-rate terms appearing in the governing 
conservation equations. 

The present calculations were also carried out to determine the 
effect of changing the initial peak frequency 6^. Figure 29 shows the 
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effect of changing the initial spatial frequency on the growth/decay 
characteristics of the coherent structure amplitude. It has been ob- 
served in non reactive calculations [35] that there exists a frequency 
for which the local amplification (-aj) is a maximum. The eigenvalue 
stability problem solved in Section 4.1 indicates the possibility of a 
peak amplification for the reactive case also. Here we present the 
solutions calculated for various initial frequencies. As the initial 
frequency is increased, the amplitude of the coherent structure also 
increases. This is again similar to the behavior observed in nonreactive 
case. The solutions for Bq = 0.4 is not shown completely except so far 
as to show the increase in |A[ . It is evident from these solutions that 
as far as the interaction is concerned, 3 q= 0.4 as compared to Bq = 0.2 
and 0.3 is one for which the initial large-scale structure is most ampli- 
fied and is thus expected that the large-scale structure at this frequency 
will be the most efficient extractor of energy from the mean flow. De- 
tailed calculations for various other Bq are however necessary to deter- 
mine which initial frequency results in peak amplification. Such calcu- 
lations have not been carried out at present. Due to the dominance of 
finite-rate kinetics, it is expected that in general solutions for various 
Bq will give very similar results as can be seen from the variation of 
the mean flow centerline specie mass fraction shown in Figure 31. There- 
fore, at present, all other reactive solutions computed here were calcu- 
lated for Bq = D.2 and numerical experiments for various Bq has not been 
attempted. 

The local eigenvalue problem, hence the local amplification rates, 
enters the interaction problem only through the integrals of the eigen- 
functions. As far as the streamwise development is concerned, because 
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of the finite amplitude of the disturbance, the local eigenvalue problem 
essentially lose their identify. Nevertheless, the local amplification 
rate, -Oj, and its streamwise development is of interest due to its ef- 
fect on the behavior of the local eigenfunctions. Figure 30 shows the 
streamwise variation of the local wave number, np and the local amplifi- 
cation rate, -aj for various initial frequency. The wave number, 
increases with streamwise distance and also increases with increase in 
frequency. This indicates that the wavelength of the disturbance 
(= l/an) decreases with increase in 6. and also decreases as the mean 
flow develops. A similar trend was observed for the non reactive case, 
also shown in Figure 30 for the case p^ = 0.3. Such decrease in scales 
were also observed by Demetriades [40] for non reactive wakes. Though 
no experimental data is available for the reactive case, it is however 
expected that the trends calculated here will also be observed in future 
experiments. The local amplification rate, -aj decreases with streamwise 
development. An interesting observation here is that on comparing the 
reactive and non reactive cases for ^^ = 0.3, the local amplification 
rate starts higher for reactions as was calculated in Section 4.1. How- 
ever, as the mean flow developes the reactive amplification rapidly de- 
creases to values below the corresponding non reactive case. There are 
two reasons for this behavior. As shown before the amplitude of the 
large-scale structure is significantly lower for the reactive case and 
also the peak value occurs around x/L = 13.2, whereas for the non reactive 
case the growth of |A|^ continues till x/L =30. This results in lower 
amplification rates ai is observed in Figure 30. With increase in 
the amplification rates also increase as is expected. However, the rate 
of decrease with streamwise development actually increases with increase 
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In initial frequency. An interpretation of these results is that in re- 
active flow, the coherent structure development is actually stabilized 
more rapidly. A similar trend has been observed experimentally by 
Ganji and Sawyer [21], Combustion therefore seems to affect the develop- 
ment of coherent structure in the mean flow in a manner consistent with 
past experiment and theoretical studies [21 -24], Essentially, the large- 
structure is initially amplified but as the interaction between the 
coherent structure and the mean flow developes, the growth rate 
(djAl^/dx) and the local amplification rate is significantly lower as 
compared to non reactive case. Detailed calculations for various other 
combustion models are necessary before any generalized conclusion can be 
drawn concerning the effect of the finite-rate kinetics on the interac- 
tion between the mean flow and large-scale structure. However, it is 
expected that results very similar to the present calculations will be 
obtained for other combustion systems. 

4.2.2 Effect of Intermittency and Molecular Diffusion 

Even in premixed turbulent flows, experimental observations indicate 
that there are situations in which the instantaneous value of fuel or 
oxidizer concentration and, therefore, the instantaneous chemical reac- 
tion rate vanish [71], Physically, this is caused by the random varia- 
tion of local flow properties coupled to high local fluctuations of the 
specie concentrations, resulting in situations in which the instantaneous 
value of either fuel or oxidizer is zero- In such situations there 
is no chemical reaction despite the finite values of average concentra- 
tions. The longer such a situation prevails, the stronger is the damping 
effect. This phenomenon has been called 'unmixedness' or 'intermittency.' 
There are several definitions of intermittency [71, 72] but since we are 
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interested in theoretical computation of finite-rate chemistry we modify 
the model developed by Spiegler et. al . [70]. This model is simple 
enough to be incorporated into the interaction problem studied here and 
is described in Appendix C. It has been pointed out that the level of 
intermittency of a certain reactant may affect the structure and behavior 
of turbulent flames and is therefore a region of great interest. The 
model described in Appendix C was incorporated into the general computer 
code and some preliminary results are presented. Figure 32 shows the 
variation of the terms contributing to the development of coherent struc- 
ture for the cases with and without intermittency. For the test condi- 
tions presented here, the effect of intermittency is to increase ttie growth 
of the amplitude. The stress production mechanism is also larger due 

to the intermittency effects. The coherent structure induced pressure 
I ^ and dissipation I ^ are also larger with intermittency. The mean 
flow dissipation I is not affected by thevinclusion of intermittency. 

The dominant production mechanism is still associated with the Reynold 
shear stress integral as one expects. Both and convert 

(the latter reversibly) disturbance kinetic energy into mean flow thermal 
energy as was observed in non reactive calculations here and elsewhere 
[39]. In the initial growth stage, the coherent structure development 
is relatively small and I primarily determines the decay rate. However, 
as the mean flow spreads and the coherent structure amplitude grows I^ 
plays a much less significant role. The shear stress production Ip^^ 
then plays a more important role in the amplitude growth rate as the 
shear layer develops. The effect of intermittency, however, is very 
marginal in the disturbance conservation mechanism essentially due to 

the fact that intermittency effects the rate constants and k. which 

i 1 
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contribute mainly to the finite-rate terms. As discussed before the 
reaction model used here is kinetically dominated and therefore not very 
significant changes in the coherent structure development is expected. 

Figure 33 shows the effect of intermittency on the other integral terms 
appearing in the conservation equations. In the linear region dissipa- 
tion I primarily determines the development of dV /dx. However as the 
coherent structure grows the diffusion flux for kinetic energy, dl|^^/dx, 
increases and is responsible for the later stages of the mean flow 
development. Also shown in this figure is the variation of total enthal- 
py diffusion flux gradient, dl^/dx. As mentioned before -dl^/dx^O, and 
the contribution of coherent structure motion to the thermal field de- 
velopment first increases and then decreases. Essentially, all coherent 
structure generated terms follow the growth/decay characteristics of 
the amplitude 1A|^. Other terms of interest are the specie diffusion 
flux integral gradient, dip, /dx which vanishes for non reactive case 

(i.e., for Y. =0). However, for reactive conditions the variation of 

these terms essentially follow the rate mechanism. Note that this term 

appearing in the specie integral equation (13) as -dl^ /dx. vJe see that 

^k 

as hydrogen and hydroxyl is consumed, +dlpv /dx and +dlp, /dx first 

oh 

increase with increase in 1 a|^ and then decreases with decay in amplitude. 

2 

dip, /dx on the other hand, first decreases and then increases with jA] . 
‘^H20 

Effect of intermittency is to increase the values of all dip /dx terms. Essen- 

k 

tially, these coherent structure induced diffusion flux gradients terms 
show contributions that increase the reduction of H 2 and OH mass fraction 
while increasing the H 2 O mass fractions. The magnitude of dip /dx is 
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temperature excess and velocity defect for the various conditions studied 
here. The effect of including diffusion is to markedly reduce the tem- 
perature increase. However, the velocity decay was not significantly 
effected by either diffusion or intermittency. The reduction in tem- 
perature excess due to diffusion effects is basically because of the 
lowered rate of change of the specie mass fraction which is shown in 
Figures 38a and b. The effect of intermittency is still almost insigni- 
ficant. The inclusion of molecular diffusion was carried out only in the 
later part of this research and due to lack of time no detailed calcula- 
tions were carried out. One major conclusion of this theoretical analy- 
sis is that the molecular diffusion term must be included for finite-rate 
calculations. 

Finally, since the results indicate that the growth rate of the co- 
herent structure amplitude during reactions is much lower than the non 
reactive case, Figure 39 presents the ratio of growth rates for the two 
cases. As expected this ratio is less than one and the growth rate ratio 
decreases with streamwise development. This is similar to the results 
in reference [21]. The increase in frequency reduced the ratio even 
further. This is attributed to the much larger growth rate for the non 
reactive case as compared to the reactive case. 

In conclusion, the interaction between the mean flow and coherent 
structure has been theoretically investigated for a supersonic non reac- 
tive and reactive wake. The non reactive solutions show excellent 
agreement with the past calculations [38, 39, 63]. For the reactive case, 
due to lack of experimental data no direct comparison was possible. 
However, the general trends of the interaction was quite similar to the 
non reactive case, except so far as to the reduced magnitude for the 
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coherent structure amplitudes. The effect of different initial ampli- 
tude, |A[p and frequency, on reactive solution has been computed. 
Turbulence was modeled in a very simple manner and included in all the 
reactive analysis. Intermittency and molecular diffusion effects were 
also included in some of the calculations. In general, the reactive 
solutions show qualitative agreement with experimental observation [21]. 

4.3 Fine-Scale Turbulence in Self Similar Mean Flow 

The turbulence model formulated in Section 4.2 was solved with the 

assumed shape functions, equation (133) for a self similar mean flow 

to determine the general trends. No detailed calculations are presented 

for this model due to lack of data necessary to evaluate the various 

constants appearing in the equations. The turbulence model utilizes a 

1 ? 

spectral partition of turbulent kinetic energy, ^ q" and turbulent to- 
2 

tal enthalpy, H" following the approach of Hanjalic et. al . [49] and 
is shown in Figure 40. In all there are 11 equations that must be 
solved simultaneously with appropriate initial conditions. The results 
presented have not been checked and their validity still remains to be 
seen. However, preliminary calculations indicate general trends similar 
to the variation observed in far field turbulence decay. This model 
was essentially formulated to account for the interaction between the 
mean flow, coherent structure and fine-scale turbulence in reactive flows 
and as such, an attempt was made to include some physical mechanism that 
is not present in classical turbulence modeling. Also, gradient approxi- 
mation was avoided except in third order correlations due to the esta- 
blished fact that in reactive flows gradient approximation has limited 
validity [3, 4]. Of course, all turbulence modeling such this one have 
some adjustable constants which have not been determined due to lack of 
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. data for wa.es. ThU .odel also incorporates a 

jetailed expenme distribution. PCO 

^ eualuate tbe effect of temperature 

,._alited temperature 

fluctuations on the mean temperature profiles and 

-^----::rr;::l--ortemperatu.^ 

probability J „ ,„d tbe mean square fluctuations. 

Ration (49) was used to ^ 

e”2 was determined by usi 

9^re 4. 9 i.es the variation , ^^duction of 

y; :,,, 1. results Show that 

hydrogen based on mean temperature, r^^ 

-■ "• '•••" '• ■•" ■•• *- 
IM ..d I. d»...l.“ .1 1 

served for other specie to 

2^ is also quite significant. 

S 43 and 44 give the Strea»ise variation of the 

Finally, Figures 43 and t,,bulent total 

, turbulent kinetic energy, Ei ) > 

centerline values reactive 

«'■' - - .... -- 

” rsntc _.*/L= 1000. 

V - 0 05 0.90, 0.025, 0.025. Xq-Xq/ 

^koo‘°’ ’ taken at present as follows. 
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. The initial conditions at x„ 

^ ‘ ' . _ / \ - in ^ 


o(Xo) = E(x„) = N(x„) = C,(x„)-10 

D„(x„) = Dm **o' ' '’s'^o' ' 
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All constants were taken to be equal to unity except for and 
= 0. 1 , k 3 * = 0.01 , and k * = = 0 

Figure 43 shows the variation of E{x)/E(x^) and N(x)/N(x^) for 
different values of correlation coefficient r^^ and Figure 44 shows 
the effect of changing (equation (71)). All solutions show decay 
characteristics as expected. The correlation coefficients r^^ are not 
constant across a developing turbulent shear layer except perhaps far 
downstream [6] and detailed numerical experiments are necessary to 
determine the effect of their variation. 



y. CONCLUDING REMARKS 


A theoretical analysis of the interaction between the mean flow and 
large-scale structure has been presented for a supersonic reactive wake. 
For numerical simplicity, von Karman integral considerations were used 
for the present solutions. The i,ntegral shape assumptions were shown to 
be quite accurate in determining the non reactive mean flow development. 
For reactive cases, the specif shape function assumptions were found to 
be quite restrictive. However, better choice of shape functions would 
depend on the availability of experimental data for the ^2~^z 
combustion system studied here. Furthermore, due to the single frequency 
model and other assumptions, point-by-point analysis of the interaction 
may not be very significant. The salient physical mechanism responsible 
for the development of the nonlinear transition region in compressible 
reactive wake has been studied for various initial conditions and the 
results are consistent with some experimental observations in other com- 
bustion studies. Improvements in the modeling of large-scale structure 
must await more refined experimental studies under well controlled con- 
conditions, especially for reactive conditions. 

Integral considerations were found to be quite restrictive on the 
shape assumptions for species but were used here for numerical simplicity. 
Finite difference calculations are, of course, possible except for the 
increase in computational costs. Since all flow fields are three dimen- 
sional, the interaction of a three dimensional coherent structure must 
be included for better understanding of the developing interaction. Fur- 
thermore, the effect of the subharmonics must also be included in future 
studies due to the observed phenomenon of pairing and breakup of large 
structure. Linear stability analysis of the reactive eigenvalue problem 
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can be extended to include both three dimensional effects and subharmonics 
but must await further guidance of more detailed experimental studies. 

The turbulence model developed here can be included into the general 
interaction problem without too many complications. However, it must 
first undergo detailed numerical testing for various combustion flow, 
which is again restricted due to lack of available data. 

The results presented here is the first theoretical attempt to in- 
clude finite-rate kinetic effects in the study of the stability and 
interaction of coherent structure in free shear layers. The inclusion 
of 'coherent structure motion in reactive flows is very important due to 
the established existence of large-scale structures in all kinds of 
shear layer and theoretical studies must therefore account for such in- 
teraction phenomenon. 
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APPENDIX A 

THE EIGENVALUE PROBLEM FOR REACTIVE FLOWS 


In this Appendix, the generalized eigenvalue problem for a two- 
dimensional reactive flow is formulated. The linearized form of the 
dimensionless equations for coherent structure motion is transformed 
using local coordinate transformation (equation (52)) and reduced to 
the stability equations by the application of the assumed shape function, 
equation (53 )• For constant transport properties the general viscous 
stability equations are 

i[(w- c)p +pu] + O)' = 0 

ia(w-c)Y|^ + pY|^'aV = Y^ ‘'■d|^ k=l,.,,v 

p(w-c)u + C^* p - i p w' V = 

ia^(w - c)v + p' = Vy ( A1 ) 


ip(w- c)[(l + V^w)V^u+ E (m,^D^T + )Y,^ 

k " 1 f 

+ M f] + i(l + V^c)C,* , + p ^H' V = Ej 

P = c* [p T 9, t f * T-y^ ^ ] 


where 


w = and -c = -11-=-^ ( A2 ) 


c 'c 

c* is the complex phase velocity with respect to the observer. The 

viscous contributions appearing on the right hand side of equations 

( A1 ) are given by the species diffusion. 


‘*k “ — [f ''t" + 2pp' Vp' - p VptpYp’ S' 
c 


'k ' 'k ■ 'k^ p 'k 

+ (p V' + 7^'):] k = 1,.^.V 
the viscous dissipations, 


( A3) 


90 




(A4 ). 


[i I v' + j av" - v] 


3 .:- 


aReV, 


(A5 ) 


and the energy dissipation 

E = [i p r + 2u' v) + p^u u" 

+ p(2pu' + ? u)u' + (7^u" + pp' u' - f ua )u] 

+ K Cj* [p^T" - T] 

V) 

C K 

D. . 

where D. , 7, F are the mean molecular diffusivity, molecular viscosity 
and thermal conductivity and are assumed known at present. Here 
prime denotes differentiation with respect to n- Furthermore, 


1 

r * s U 

2 

V 


and 


C^* = l/(Y-l)Me^ReLPr 


( A6) 


and all other terms are defined in the text. 
i) Inviscid Solution 

The inviscid equations are obtained from equations (Al ) by 
neglecting all terms containing diffusivity, viscosity and thermal con- 
ductivity. The resulting equations can be manipulated to obtain a second 

order differential equation in terms of the pressure perturbation p(n) • 

2 2 

R(n) = 0 (A7) 


where 


R(n) = aV T 2(w - c) I ^ ^ 

^ k=l MT T 

Here G. corresponds to the kth specie production rate due to instability 

lx 
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wave. It can be shown that for nonreactive flow equation (A7) reduces to 
the classical pressure perturbation equation of Lees and L1n [43]. How- 
ever. for reactive flows the tern. Is determined such that 

G, ^31*^ .^41^‘^ ^5i3 

G2= T,2t»2 ^22*^3 ^32*5 ^52^ ^ , 

G3=p'['v, T,3.Y2 T23ti3 TjjtY, T„tp T„] 

G, = p^[il 1,4*^3 ^34*'^ ^44*‘> V 

where G.(i = l,4) and V, (1 = 1.4) indicate species Hj. O 2 . OH and H 2 ) 

^ -TM-itii-lA) are coefficients which are 

respectively. The matrix T,j(i - 1.5 J 1.4) ore 

r AQ^ in ( A8) and on some algebraic 
described below. On using equation ( A9 ) in 1 «o , 

manipulations we get ^ ^ _ 

p"vS’[-|?^tn„.-i)-V-6] 


(AlO) 


v/(w-c) V* 1 1 


(^.g)-l]oVT2 = 0 


where v* and I are complex constants appearing due to finite rate kine- 
tics and are defined as 


= [G^. i i v/ 


E(U - c)aT 

i2 


'I 


-2 « 

and 2 ,- \2 ^ *^o T TT 4 . ^ 7' 

. , [G, * 1 "r ■* ’ ®l 

T|?V E(u-c)aT 

where prime denotes differentiation with respect to n and 


(All) 

(A12) 
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6o - ’ H j, V "k 
'=1 ■ I ^k* «k 


(A13) 


6 „ - Z N( , Nn - 1 

0 j,, 0 J 6 , 


(A 14 ) 


and 


Bi * S fiL M. 
1 k=i 


E ■ , "'k “k \ 

k= 1 


(A15) 


The constants C^, C|^* and Dj^ are defined in the text. The terms in (A13) 
and (A14) are due to the algebraic manipulation involved and are given by 


k^l’.A 


(A16) 


The matrix T^j^ which also appears in equation (A9) is the finite rate pro- 
duction of each specie due to the disturbance field and is given as 


T,2 ■ *2 * \ 


Tj2 - T<2 a T, 


^32 ■ - 2 ’'2 ® \ ^3 


T42 ■ 0 


^52 ' ^ *2 ^ ^2 ' ® '‘ b , 


93 



^14 ■ *4 C P “fj ^ 


T 24 -O 


T34 - 2 X4 c r T, Yj 


^44 ■ - 2 *4 B kt,^ T 3 

T 54 ■ ><4 [3 C k,^ Y, Y^ - i 0 Y|] 

T * X r-^ + 

'll 


T21 « T22/X2 


'*^31 “ ^1 

^41 “ 1 ^ 44 /^ 

^51 ' " ^ 54 ^V 

T^3 » X3 [T^2/X2 - T^^/X^l 

^23 " ^3 ^^22^^2^ 

"^33 ' ^3 ^"^32^^2 " ”’^34^^4^ 


^43 ' “ ^3 "44^^4 


T53 = X3 [T52/X2 - T54/X4] 


(A17) 
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where 


and 


C « P. 


X, « 


'^2e 


M2 

f 

!k 

Y 2 
^3e 

' «, 

M2^ 

Ae*- 

Ml 

“e 


p^L 

Mo 


L 

s 

Y2e 

Pc'- 

M3 

“e 

Y3e 

Pe*- 

M4 

“e 



B 


M, 


4e 


(A18) 


(A19) 


where are the 1-th species molecular weight and subscript e indicates 
reference edge values. 

The terms N|^ and appearing in (A13) are further defined as follows, 

2 1 2 

Ni = + N2 


N« » N 


1 


'1 


' ^1^ ^1 ^3^ 

N4 ' + N2^ N2 + N 3 ^ 
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1 2 2 

• M 

M2 - M^ 

M3 « M^ + N2^ M2 + M^ 

M4 “ M^ + N2^ M2 + M^ 

[T ,2 ^2 ■" ^Z2 '^52^'''^ 

D = r(l-'<'i T22 )->^i T^2 ^1^-’''1^32 
M - T ^2 ^^■*■''’32 " ’*'2 

N,2 . ♦, ITj, + T 3 , + T,, 

' *1 C^51 * ^31 * ^41 **3 ^^'’2 

M^ ' ['*', tT3^ M^ + M^} - tl'2 ^i' 

■>2 ' (’ - +, T„) - 4, Hj, n/ 4 T,, n/) 

N,3 - ♦, IT , 3 4 T„ N,'>]/(1 - ♦, T 33 ) 

- ♦, IT 23 + T„ N^'^/O - ♦, T 33 ) 

' *1 ^^53 * ^43 ^3 ' *1 ^33^ 

■ I- *2 T3' P' 4 *, T„ M^]/(l - ♦, T33) 


(A20) 


(A21) 


(A22) 


(A23) 
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^ 1 ^ “ ‘ * *1 ^34 ^ 13 ^ 

^^2^ • \ "^34 "^32 

^ 3 ^ * dJ ’ ■ ’^1 ^ 33 ^ ^54 *1 ^34 ^ 53 ^ 

A 

to - ♦, T33) Y^' * ♦, T3, 73'] 

4 

D^ - [(1 - »P-| '>'33) n - *1 T^4) - ''’34 ''^43^ 


where 


’*'1 ‘ i a(w-cjv: 


C-j p 


'^Z * — \2»T 

^ a (w- c) V. 


(A25) 


ii) Viscous Solution 

The full viscous equations for this stability problem equation 
(Al) is very complicated and no solution has been obtained so 
far. Lees and Lin [43], Dunn-Lin [66] and Lees and Reshotko [58] 
have in the past presented solutions for the viscous equations by solv- 
ing a set of reduced equations that retain terms up to a certain order, 
either near the critical point or near the surface. Dunn-Lin [66] 
ordering of the various terms accounts for supersonic flows and for the 
propagation wave velocity c to be a substantial portion of the free- 
stream velocity. Therefore, at present a similar ordering is applied 
here whereby 
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(j), 4 ,' % ©( 1 ) 

^ o(^) , (W- c) -v. 0(1) 

u 'V- 0(1) , V 'V- 0(eu) ,T'\-u .Yj^-v-u ,p-vu (A26) 

p 'V u 

^ ' “k 

Using these ordering, equation (Al) can be reduced to a form similar 


to Dunn-Lin viscous equations 

1 Vp _ - 

u"' - i a Re V ‘ (w- c)u' = 0 

p ^ 

;. + i ^ + j I 

p p T k=l 

aV 1 Yk* 

^k" " ~2 ' 'TJT^ ’ k=i,...v 

k p ^ p 

T" - i(-^) (w- c)MT 

K Cx* P 


(A27) 

(A28) 

(A29) 


= (• 


V 

I 


1 


K Cj* 


k=l a p 




(A30) 


+ (i(^)^(w-c) I (m,^D,^T+E,^ -h|^)Y,^ 

kCj p k=l T 

We see that this is a eighth-order system of ordinary differential 
equations dependent basically on one parameter (a Re). For nonreactive 
case, the equation (A29) and all terms containing Y|^ drop out and the 
resulting equations are then identically Dunn-Lin [66] viscous equa- 
tions. Lees-Reshotko [58] ordered terms in two parameters (otRe) and 

M ^ and retained additional terms containing temperature and viscosity 
e 

fluctuations. Such an ordering would eventually be necessary to deal 
with viscous reactive stability but have been presently neglected by 
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assuming all transport properties are function of mean variables only. 

Of course, equations (A27) to (A30) can be directly solved simultaneously 
with appropriate boundary conditions. However, due to the analytical 
nature of these equations, some modification is possible. Consider the 
behavior of these equations in the outer flow where all mean flow varia- 
bles have attained their external values. A solution for equation (A27) 
in the outer flow can be immediately written as 

u' (n) = exp[± i{i aRev”^ V^cl'^n] . n -► (A31 ) 

However all other equations are coupled and no direct solution can be 
written for the general reactive case. This is in contrast to nonreac- 
tive case forY|^* = 0) when we get a solution for equation (A28) as 

Y|^ 'V' exp[± i{i a D|^ ^ cl'^n] > n ->■ (A32) 

and when Yj^ terms are neglected, equation (A30) gives 

T -x. exp[± i{iCic Cj*)"^ac V^M}\] , n -*■ ±® (A33) 

Solutions (A31) and (A33) are identical to the ones obtained by 
Lees-Reshotko [58]. These are all linearly independent solutions and 
governing equations must satisfy these solutions in the outer edge. How- 
ever, for reactive wakes only equation (A31) is a direct solution since 
all other equations are coupled. For n ^ the positive exponent cannot 
satisfy the outer inviscid boundary conditions and must be neglected. 

Then on applying transformations [58] 

h=v/u , f=u7u , e = T/u and \'Y,^/u , k = l,...v (A34) 


we get 


h' 

f" 


s f - - - ^ 6 + Z 

p FT k 

„ i a Re V ^ 

= -f(r + 3f')+ — T 2 — V^(w-c)f 

p 


(A35) 

(A36) 
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II 


iaV, 


L " = “2f ii|^' - (t +f )£|^ + ;— . k = 1,...V 

p 

e"«-2fe'- (f +r)0 + [iav^(w-c){ Z (m,^DJ (A37) 

(p7C^*) k=l 


* ’’kSk^ 


An approximate set of boundary conditions at 
T^-^oo f^T^-k-oo) = -i[i a Re V ^ V c]^ 




1 


[p{iaRe V ^ V^.c}^] 


* 0 

6(it^ “) = 0 


For anti syimietri cal disturbance 
f(0) = h(0) = 0 
e' (0) = £|^' (0) = 0 


(A38) 


(A39) 


For the general solution of viscous reactive flow stability, the 
equations (A35) to (A37) must be solved for the given boundary conditions 
such that the outer solution satisfies the inviscid conditions. Essen- 
tially two characteristic values a and Re must be determined for the gen- 
eral viscous case in contrast to only a for the inviscid part. Neutral 
stability characteristics can then be obtained as a function of (a, Re) 
for a given Mach number. This generalized stability analysis of viscous 
reactive wakes will be an interesting research area for future study. 
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APPENDIX B 


INTERACTION INTEGRALS FOR COHERENT STRUCTURE MOTION 


The Integrals appearing in the governing conservation equations in 
Chapter II can be rewritten using the assumed mean flow profiles and 
using the coherent structure and fine-scale turbulent representations. 
Since at present only the interaction between mean motion and coherent 
structure has been studied therefore only the Integrals relevant for 
this analysis Is presented here. 

Species diffusion flux: 


= 26lAl^[V - l)Re(pG*)+Re(OY *)}dn 

^ -00 ^ 

00 

+ / ■^Re(Y|^p*)dn] t k=l,...v 


(B1 ) 


Species production: 

OD '■ - - 

Ru = 6 / r ° dn 


k = 1 , . • . V 


(B2 ) 


where rj^° is the total mean production rate for the k-th species and can 
be rewritten using equation (140) as follows 


t Y„„21;|2)) - (A/o/2|A|^(p^lie(Yo^*) 

* 2 7 Y„^Re( * 2 Re(; 9„^*) ^ I'p I ^) )] 


(B3 ) 




* ''hj ^ ''oh''h2 ""‘o ''oh*) * ''oh‘ ■> 

1 2 ,—2 1 G /M 2 


T7 2 


* 6 | Al 2 pY „^ Y (,^ 2 | g | 2 ,. 0 Fj ,^( p " Y „^0 + 2 | A |‘( p ‘| Y „^ 0 r 

,2 


7 

w - 


(B4 ) 




101 



0 


= E r ° 1 F r ° 




( B5) 
( B6) 


where all the constants are defined in the text (Section 2.7 ) and * de- 
notes the complex conjugate, and k|^ are the mean rate constants 

which for this interaction problem is the function of the mean tempera- 
ture. Here Re stands for the real part of the complex product. 


Mean flow total enthalpy flux: 

1^=1 ^ + I ^+I ^ 

E 'Ei '^2 ^3 


where 


^ = 26lA|^V^ / i (H- 1) Re (p u *)dn 

1 -oo P 

^ = 26|Al\ / [V^u lul^ + MRe(uf *) 

2 -00 

+ z(m^DJ+e,^ )Re(uY^*)]dn 

k ^ 

oo 

Ip = 26lAl^ / [V^uRe(u p*) + MRe(pf *) 


+ z(m^D^T+e,^ )Re(p Y,^*)] z 
k f p 


( B7 ) 
(B8) 


(B9 ) 


(BIO) 


k T p 

The integral for disturbance generated diffusion flux in the mean flow 
kinetic energy conservation is 


^ 2 

^kE^ " ^ Re(pu*) + 2V^ulul^]dn 


(BID 


The Reynolds stress production In ^ is given by 

Ino^ = Inn + In. + I n. + I nn ^ + I , ^ 


RS 


'RS, 


‘RS. 


‘RS, 


‘RS, 


‘RS^- 


(B12) 


where 
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= - 2 lA]^V^^ / pRe(ua* V*) dn 


1.5 ^ = -26|Ai\ / ^Re(pu*)(|^) dn 
^2 -® p y 

IrS,' ‘ -2«l''l V / 


'3 

c 


'rs/ ‘ / vRe(pu*) l^dn 


I,s; = -2|AiV / Pl-l'|^dn 


where 


f 9 \ _ / 9 \ _ JT_ d6 9 
^9X\. ~ ^9X' ~ 6 ^ 9n 


The pressure work integrals are 

d 


Ip =-C, ^ 6 I (p-l)^dn 


'1 


I =C, 6 ; (p-i,[(B) .£|ljdn 

2 y 

The mean flow pressure is determined in general, from 

P = 1 - V [271A1^ |av|2] 

^1 

The organized motion pressure work integral is 

_ 00 ^ ^ 

I^ = -2C-|6 |Al V / Red apu* + avp^*}dn 

P * ^ -00 

The mean flow viscous dissipation integral is 

I ' f — 

° 5^ i " " ” 

where, for a linear viscosity- temperature relation, U=CT, C = 1 
p y = 1 - 2 I A I ^[Re(p T*) + ^ |a V | ^] 

The organized motion viscous dissipation integral is 


(B13) 

(B14j 

(B15) 

(B16) 

(B17) 

(B18) 

(B19) 
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(B20) 


2v/|Al^C - 


f / t 


2 ^3 


iur + 


+ (4 |av’ 1^ + lu' 1^) +1 {Re(2 i V u'- 


4 i jciil^ uv'*)} ]dn 


The relationship between shear layer thickness and velocity defect is 
given by equation ( 48 ) which results in 


6 = 


Cd/2 


Vc(Ci-VcC,) 


(B21) 


where 81 = /^ and = >^^/2 and Cp is the drag coefficient of the body. 


The disturbance kinetic energy integral is 


6 

7 


/ u( u'^ +v-^)dn = 6 |A 1 ^I^ -V^I^ ) 

.00 12 


(B22) 


where 


= / ( I u I ^ + jciv j^)dn = 1 

1 -00 

Ik = /(|u|^|<.v|2)e-''^<‘" 

2 


(B23) 


The filtering procedure: 

On using equation (10) the mean value is given by 

i'lli = (ij; + (f' + + I*'* + + r'<p' (B24) 

On subtracting the mean value from the total value and taking condi- 
tional average we get 

<44 - + (f-V + R (B25) 

V V 

Fine-scale turbulence contribution is obtained by 

<}>ip - - (|4 = (B26) 

Similar filtering is also applicable to the 0(3) or higher fluctua- 
tion correlation. 
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APPENDIX C 

MODELS FOR TURBULENCE. INTERMITTENCY AND DIFFUSION 
Turbulence Models 

The algebraic eddy viscosity models used in these calculations were 
modified from the formulation of Spalding [69] and Schlichting and Hinze 
[70]. They are given in dimensional form as follows. 

Model I [69] 

„t* = 0.037 6*1 (Cl) 

where U„,„ and U„- are the local maximum and minimum values of the stream- 
max min 

wise velocity and is the dimensional turbulent viscosity. 

Model II [70] 

p^* = 0.044 a* (C2) 

where the characteristic density thickness a* is given in the dimensional 
form as 

A* =1 /F* (u_^- U*) dy*| 

— oo 

On nondimensionalizing equation (Cl) we get 

€j = 0.037 6 p (C4) 

where 4 T = U 4 .*/p U L is then idenfitied as the dimensionless turbulent 

1 t oo CO 

eddy viscosity. Equation (C2) then becomes 

Cjj = 0.044 a (C5) 

where the nondimensional density thickness 

A = I / F(1 -F)dy| = /if 6 (C6) 

— oo 

and ^TT *” 1 ^ 4 - /p ^ 

11 t 0° ® 

During the calculations it was found that the turbulent eddy 
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viscosity Model II gave slightly better results and was used for all the 
reactive calculations. 


Intermittency Models 

The intermittency model used here is basically the model proposed 
by Spiegler et. al, [71] and is repeated here for completeness. This 
model is applied to the reaction 


H 


2 



2H2O 


(C7) 


and the determination of the intermittency factors for the forward and 
backward rates are given by the following steps. 

For each specie (k), the mean square fluctuation of mass fraction 

— 2 

Y|^ is determined from the dimensional relations 


7T 


= Cj 6*" 


A 

9y* ' 3y* 


20.25 


(C8) 


We then compute 



> = 

^k y 2 
'k 


(C9) 


and 

A. - 0.233 

Then we determine the intermittency factor for reaction (i) as follows: 

— 2 

U^. =max(Uj, ), if the respective Y|^ 's have the same sign. 

~ k^i respective Y^'s have different signs. 

The finite-rate constants for each reaction (i) and direction 
(m=f,b) then becomes 


k 


mi 



(1-U.) 


(Cll) 
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These rate constants are then used in the finite-rate production terms 
appearing in the governing equations. 

Diffusion Model 

The diffusion term due to molecular diffusivity in the integral 
form is given in the nondimensional form as 

" _ - 2 “577-2 

'd^* ' J » (C12) 

Applying the transformations we get 

= } / P^D^[Y|^'^+ 2|A|^|Y^'|^]dn (C13) 

where prime denotes differentiation with respect to n. For simplicity, 

we relate the mean molecular diffusivity to the kinematic viscosity by 

the Schmidt Number Sc (= v/D|^) which at present we take to be unity. 

This term is then added to the right hand side of equation (13) 
k 

which results in an additional term in the thermal energy balance by 
virtue of the heat flux term equation (8). 
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igure 1 The geometry of the physical problem: Reacting free shear layer 
in the wake of a flat plate. 


Figure 2 
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The path of numerical integration in the complex 
n-plane for the reactive eigenvalue problem. 
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Figure 3 Neutral stability characteristics for non-reacting 
flow: critical point, wave number and relative 
phase velocity as a function of temperature excess 
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Figure 6 Spatial neutral wave frequency as a function of relative Mach number 
and wake temperature excess. 
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Figure 7 The eigenvalues for symmetrical oscillations as a function of wave frequency. 
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Figure 9 The eigenvalues for symmetrical oscillations as a 
function of wake velocity defect. 
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Figure 10 The eigenvalues for symmetrical oscillations as a 
function of wake temperature excess. 
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Figure 12 The eigenvalues for anti -symmetrical oscillations 
as a function of free stream Mach number. 
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Figure 13 The eigenvalues for anti -symmetrical oscillations as 
a function of wake velocity defect. 
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Figure 15 The characteristic variation of the pressure eigen- 
function amplitude [pj and phase <j>p across the 

shear layer. For symmetrical oscillation, 
p'(0) = 0 and M„=2, g = 0.1, T = 1500°K. 
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Figure 16 The characteristic variation of the pressure eigen- 
function amplitude |pl and phase 4)p across the shear 

layer. For anti -symmetrical oscillation, p(0) = 0 
and M = 2, 6 = 0.1, T = 1500°l(. 
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Figure 17 The characteristic variation of the amplitudes of pres- 
sure eigenfunction lp(n)| and vertical velocity eigen- 
function |v(n)l across’ a symmetric jet. For symmetri- 
cal oscillations, M =2.0, T = 1500°K, T =3.0 and 
B=0.1. “ " 


130 



Figure 13 Phase changes of the pressure eigenfunction ;j> and 

P 

vertical velocity eigenfunction (f>^ across a symmetric 
jet. Test conditions same as in Figure 17. 
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(BATT AND KUBOTA [?7]) 
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(a) Hydrogen mass fraction profiles, Y„ (b) Oxygen mass fraction profiles, 
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Figure 21c Hydroxyl mass fraction profiles. 
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TURBULENCE: ALGEBRAIC 
eddy VISCOSITY MODEL II 
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Figure 24 Development of the average coherent 
structure energy density, C, stress 
dissipation I^^ and mean flov; 
dissipation I^ in the turbulent 
supersonic non reactive wake. Flat 
olate at M =3.0, Re = 70 ,000/cm, 

» to 

g_ = 0.3. Cn = 0.031, x^=5.0, = 0.2 






Figure 25 


Gro'v/th of centerline mean flov/ temperature excess for 
various initial amplitudes of disturbance. Turbulent 
reactive test conditions: M^=3.0, Re = 70 ,000/cm, 

T^= 1500°K, L =0.3cm, 3^=0. 2, Xq = 5.0, =0.23 

= 1.275 , (Y. = 0.0, k = 1,4). ° 






defects hydrogen and hydroxyl nass fraction 


(b) The centerline mean flow water vapor mass fraction excess. 
Figure 21 of the centerline mean flow specie mass 
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Growth and decay characteri sti cs of the amolitude of coherent structure for 
various initial values. Turbulent reactive test conditions as in Figure 25. 
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f^igure 32 


Development of the coherent st 
mecnanism •.;lth and -.vithout Tnt, 
reactive test conditions as in 


structure energy 
ntermittency. ' 
Figure 25. 
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Figure 33 Development of the centerline mean flow specie mass frac- 
tions in turbulent reactive wakes. Effects of intermittency 

o 

and diffusion for |Al ‘‘ = 4x10 . Test conditions as in 

Figure 25. 

(a) The hydrogen and hydroxyl mass fraction defects. 

(o) The water vapor mass fraction excess. 
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